function transect(ICQ4,fem_grid_struct,xt,yt,dinter,sfield,vecflag)
%TRANSECT - compute and draw transect from a fem_icq4_struct
% 
%  Inputs: fem_icq4_struct - .icq4 data structure (READ_ICQ4,ICQ4INFO),
%                            as in "icq4=read_icq4('ITER_03820.icq4',1);"
%          fem_grid_struct - FEM grids structure (LOADGRID)
%          xt              - x-coordinates of transect ([x1 x2]) 
%          yt              - y-coordinates of transect ([y1 y2])
%          dinter          - distance interval to interpolate over         
%          sfield          - the scalar field to contour.  Must be
%                            a string in single quotes, one of the 
%                            following:
%                             'UZMID'  
%                             'VZMID'
%                             'WZMID'
%                             'Q2MID'
%                             'Q2LMID'
%                             'TMPMID'
%                             'SALMID'
%                             'u'  for the normal flow component
%          vecflag         - 0|1 flag to plot parallel velocity vectors
%
%  Call as: transect(ICQ4,fem_grid_struct,[x1 x2],[y1 y2],dinter,sfield,vecflag)
%  Example: transect(ICQ4,b150,[x1 x2],[y1 y2],2000,'TMPMID',0)
%
%
% Outputs: Color surface drawn to current axes
% Little error checking is done.
%

if nargin~=7
   disp('Wrong number of arguments to TRANSECT.  Must have 7!!')
   disp('transect(ICQ4,fem_grid_struct,[x1 x2],[y1 y2],dinter,sfield,vecflag)')
   return
end

if ~(strcmp(sfield,'UZMID')|strcmp(sfield,'VZMID')|strcmp(sfield,'WZMID')|...
     strcmp(sfield,'Q2MID')|strcmp(sfield,'Q2LMID')|strcmp(sfield,'TMPMID')|...
     strcmp(sfield,'SALMID')|strcmp(sfield,'u'))
   error(['field name of ' sfield ' not valid.'])
end

% Check dims of [x1 x2],[y1 y2];
[m,n]=size(xt);
if m~=1 | n~=2
   error('[x1 x2] not 1X2')
end
[m,n]=size(yt);
if m~=1 | n~=2
   error('[y1 y2] not 1X2')
end

R=searibbon(ICQ4,fem_grid_struct,xt,yt,dinter);
ind=find(~isnan(R.UZMID(:,1)));

%%%%%%%%%%%%%%  ROTATE VELOCITY INTO TRANSECT PLANE %%%%%%%%%%%%%%%%
tgv=atan(diff(yt)/diff(xt));
u=R.UZMID(ind,:);
v=R.VZMID(ind,:);
[u,v]=rot(u,v,tgv);
R.u=u;
R.v=v;

% U POSITIVE TO RIGHT OF TRANSECT DIRECTION (OUT OF SCREEN)
% V POSITIVE ALONG TRANSECT

SCALAR=reshape(getfield(R,sfield),size(R.ZMID(ind,:)));

pcolor(R.d(ind,:),R.ZMID(ind,:),SCALAR)
xlabel('Distance along transect (meters)')
ylabel('Depth (meters)')

colormap(jet(16))
shading interp

curax=gca;
hcb=colorbar;
axes(hcb)
%ylabel('Through-Transect Velocity (m/sec)')
ylabel('Temperature')
axes(curax)

if vecflag
   vec_scale=.1;
   ind1=find(~isnan(R.UZMID(:)));
   hvec=vecplotstick(R.d(ind1),R.ZMID(ind1),v,R.WZMID(ind1),vec_scale,'k',...
      10,'m/s');
end
   
line(R.d(:,1),R.ZMID(:,1),'Color','k','LineStyle','-','linewidth',2)
drawnow

