# VMT_GridData2MeanXS

## PURPOSE

Generates a uniformly spaced grid for the mean cross section and

## SYNOPSIS

function [A,V] = VMT_GridData2MeanXS(z,A,V)

## DESCRIPTION

     Generates a uniformly spaced grid for the mean cross section and
     maps (interpolates) individual transects to this grid.

     (adapted from code by J. Czuba)

     P.R. Jackson, USGS, 12-9-08

## CROSS-REFERENCE INFORMATION

This function calls:

This function is called by:

- VMT_ProcessTransectsV3_new Not used by current version
- VMT_ProcessTransectsV3_new_display Not used by current version

## SOURCE CODE

    0001 function [A,V] = VMT_GridData2MeanXS(z,A,V)
    0002 % Generates a uniformly spaced grid for the mean cross section and
    0003 % maps (interpolates) individual transects to this grid.
    0004 %
    0005 % (adapted from code by J. Czuba)
    0006 %
    0007 % P.R. Jackson, USGS, 12-9-08
    0008
    0009 %% User Input
    0010
    0011 xgdspc = A(1).hgns; %Horizontal Grid node spacing in meters  (vertical grid spacing is set by bins)
    0012
    0013 %% Determine uniform mean c-s grid for vector interpolating
    0014 % Determine the end points of the mean cross-section line
    0015 % Initialize variable with mid range value
    0016 V.xe = mean(A(1).Comp.xm);
    0017 V.ys = mean(A(1).Comp.ym);
    0018 V.xw = mean(A(1).Comp.xm);
    0019 V.yn = mean(A(1).Comp.ym);
    0020
    0021 for zi = 1 : z
    0022
    0023     V.xe = max(max(A(zi).Comp.xm),V.xe);
    0024     V.ys = min(min(A(zi).Comp.ym),V.ys);
    0025     V.xw = min(min(A(zi).Comp.xm),V.xw);
    0026     V.yn = max(max(A(zi).Comp.ym),V.yn);
    0027
    0028 end
    0029
    0030 % Determine the distance between the mean cross-section endpoints
    0031 V.dx = V.xe-V.xw;
    0032 V.dy = V.yn-V.ys;
    0033
    0034 V.dl = sqrt(V.dx^2+V.dy^2);
    0035
    0036 % Determine mean cross-section velocity vector grid
    0037 V.mcsDist = linspace(0,V.dl,floor(V.dl/xgdspc));                                  %%linspace(0,V.dl,V.dl); Changed to make it user selectable (PRJ, 12-12-08)
    0038 V.mcsDepth = A(1).Wat.binDepth(:,1);
    0039 [V.mcsDist V.mcsDepth] = meshgrid(V.mcsDist,V.mcsDepth);
    0040
    0041 % Determine the angle the mean cross-section makes with the
    0042 % x-axis (longitude)
    0043 % Plot mean cross-section line
    0044 if V.m >= 0
    0045     V.theta = atand(V.dy./V.dx);
    0046
    0047     figure(1); hold on
    0048     plot([V.xe V.xw],[V.yn V.ys],'ks'); hold on
    0049
    0050     if V.mfd >= 270 | V.mfd < 90 %Flow to the north
    0051         V.mcsX = V.xw+V.mcsDist(1,:).*cosd(V.theta);            %
    0052         V.mcsY = V.ys+V.mcsDist(1,:).*sind(V.theta);
    0053
    0054     elseif V.mfd >= 90 & V.mfd < 270 %Flow to the south
    0055         V.mcsX = V.xe-V.mcsDist(1,:).*cosd(V.theta);            %
    0056         V.mcsY = V.yn-V.mcsDist(1,:).*sind(V.theta);
    0057     end%
    0058
    0059     plot(V.mcsX,V.mcsY,'k+'); hold on
    0060     plot(V.mcsX(1),V.mcsY(1),'y*'); hold on
    0061
    0062 elseif V.m < 0
    0063     V.theta = -atand(V.dy./V.dx);
    0064
    0065     figure(1); hold on
    0066     plot([V.xe V.xw],[V.ys V.yn],'ks'); hold on
    0067
    0068     if V.mfd >= 270 | V.mfd < 90 %Flow to the north
    0069         V.mcsX = V.xw+V.mcsDist(1,:).*cosd(V.theta);            %
    0070         V.mcsY = V.yn+V.mcsDist(1,:).*sind(V.theta);
    0071
    0072     elseif V.mfd >= 90 & V.mfd < 270 %Flow to the south
    0073         V.mcsX = V.xe-V.mcsDist(1,:).*cosd(V.theta);
    0074         V.mcsY = V.ys-V.mcsDist(1,:).*sind(V.theta);
    0075     end%
    0076
    0077     plot(V.mcsX,V.mcsY,'k+'); hold on
    0078     plot(V.mcsX(1),V.mcsY(1),'y*'); hold on
    0079
    0080 end
    0081
    0082 V.mcsX = meshgrid(V.mcsX,V.mcsDepth(:,1));
    0083 V.mcsY = meshgrid(V.mcsY,V.mcsDepth(:,1));
    0084
    0085 clear zi
    0086
    0087 %% Determine location of mapped ensemble points for interpolating
    0088 % For all transects
    0089
    0090 %A = VMT_DxDyfromLB(z,A,V); %Computes dx and dy
    0091
    0092 for zi = 1 : z
    0093
    0094     % Determine if the mean cross-section line trends NW-SE or SW-NE
    0095     % Determine the distance in radians from the left bank mean
    0096     % cross-section point to the mapped ensemble point for an individual
    0097     % transect
    0098    A(zi).Comp.dx = abs(V.xe-A(zi).Comp.xm);  %This assumes the easternmost bank is the left bank--changed PRJ 1-21-09 (now uses VMT_DxDyfromLB--not working 2/1/09)
    0099
    0100     if V.m > 0
    0101         A(zi).Comp.dy = abs(V.yn-A(zi).Comp.ym);
    0102
    0103     elseif V.m < 0
    0104             A(zi).Comp.dy = abs(V.ys-A(zi).Comp.ym);
    0105
    0106     end
    0107
    0108     % Determine the distance in meters from the left bank mean
    0109     % cross-section point to the mapped ensemble point for an individual
    0110     % transect
    0111     A(zi).Comp.dl = sqrt(A(zi).Comp.dx.^2+A(zi).Comp.dy.^2);
    0112
    0113     % Sort vectors by dl
    0114     A(zi).Comp.dlsort = sort(A(zi).Comp.dl,'ascend');
    0115
    0116     % Map indices
    0117     for i = 1 : A(zi).Sup.noe
    0118         for k = 1 : A(zi).Sup.noe
    0119
    0120             if A(zi).Comp.dlsort(i,1) == A(zi).Comp.dl(k,1)
    0121                 A(zi).Comp.vecmap(i,1) = k;
    0122
    0123             end
    0124         end
    0125     end
    0126
    0127     % GPS position fix
    0128     % if distances from the left bank are the same for two ensembles the
    0129     % the position of the right most ensemble is interpolated from adjacent
    0130     % ensembles
    0131     % check for repeat values of distance
    0132     sbt(:,1)=diff(A(zi).Comp.dlsort);
    0133     chk(1,1)=1;
    0134     chk(2:A(zi).Sup.noe,1)=sbt(1:end,1);
    0135
    0136     % identify repeat values
    0137     A(zi).Comp.sd = (chk==0) > 0;
    0138
    0139     % if repeat values exist interpolate distances from adjacent ensembles
    0140     if sum(A(zi).Comp.sd) > 0
    0141
    0142         % bracket repeat sections
    0143         [I,J] = ind2sub(size(A(zi).Comp.sd),find(A(zi).Comp.sd==1));
    0144         df=diff(I);
    0145         nbrk=sum(df>1)+1;
    0146         [I2,J2] = ind2sub(size(df),find(df>1));
    0147
    0148         bg(1)=(I(1)-1);
    0149
    0150         for n = 2 : nbrk
    0151             bg(n)=(I(I2(n-1)+1)-1);
    0152         end
    0153
    0154         for n = 1 : nbrk -1
    0155             ed(n)=(I(I2(n))+1);
    0156         end
    0157
    0158         ed(nbrk)=I(end)+1;
    0159
    0160         % interpolate repeat values
    0161         A(zi).Comp.dlsortgpsfix = A(zi).Comp.dlsort;
    0162
    0163         for i = 1 : nbrk
    0164             for j = bg(i)+1 : ed(i)-1
    0165                 % interpolate
    0166                 if bg(i) > 0 && ed(i) < length(A(zi).Nav.lat_deg)
    0167
    0168                     den=(ed(i)-bg(i));
    0169                     num2=j-bg(i);
    0170                     num1=ed(i)-j;
    0171
    0172                     A(zi).Comp.dlsortgpsfix(j,1)=...
    0173                         (num1/den).*A(zi).Comp.dlsort(bg(i))+...
    0174                         (num2/den).*A(zi).Comp.dlsort(ed(i));
    0175
    0176                 end
    0177
    0178                 % extrapolate end
    0179                 if ed(i) > length(A(zi).Nav.lat_deg)
    0180
    0181                     numex=ed(i)-length(A(zi).Nav.lat_deg);
    0182
    0183                     A(zi).Comp.dlsortgpsfix(j,1)=numex.*...
    0184                         (A(zi).Comp.dlsort(bg(i))-...
    0185                         A(zi).Comp.dlsort(bg(i)-1))+...
    0186                         A(zi).Comp.dlsort(bg(i));
    0187
    0188                 end
    0189             end
    0190         end
    0191
    0192     else
    0193
    0194         A(zi).Comp.dlsortgpsfix = A(zi).Comp.dlsort;
    0195
    0196     end
    0197
    0198     % Determine velocity vector grid for individual transects
    0199     [A(zi).Comp.itDist A(zi).Comp.itDepth] = ...
    0200         meshgrid(A(zi).Comp.dlsortgpsfix,A(zi).Wat.binDepth(:,1));
    0201
    0202     clear I I2 J J2 bg chk df ed i j nbrk sbt xUTM yUTM n zi...
    0203         den num2 num1 numex
    0204
    0205 end
    0206
    0207 clear zi i k check
    0208 %% Average ensembles from individual transects using a spatial average to points on the uniform mean c-s grid
    0209 % Fill in uniform grid by averaging ensembles (of individual transects mapped onto the mean
    0210 % cross-section) within one-half a grid spacing on either side of a uniform mean
    0211 % c-s node. This uses all of the ensembles.
    0212
    0213 for zi = 1 : z
    0214 % reorder the ensembles so the distance increments across the c-s
    0215 A(zi).Clean.bsvecmap = A(zi).Clean.bs(:,A(zi).Comp.vecmap);
    0216 A(zi).Clean.vDirvecmap = A(zi).Clean.vDir(:,A(zi).Comp.vecmap);
    0217 A(zi).Clean.vMagvecmap = A(zi).Clean.vMag(:,A(zi).Comp.vecmap);
    0218 A(zi).Clean.vEastvecmap = A(zi).Clean.vEast(:,A(zi).Comp.vecmap);
    0219 A(zi).Clean.vNorthvecmap = A(zi).Clean.vNorth(:,A(zi).Comp.vecmap);
    0220 A(zi).Clean.vVertvecmap = A(zi).Clean.vVert(:,A(zi).Comp.vecmap);
    0221 A(zi).Clean.depthvecmap = nanmean(A(zi).Nav.depth(A(zi).Comp.vecmap,:),2)';
    0222
    0223 % determine one-half grid spacing on either side of a node
    0224 plus=V.mcsDist(1,:)+(V.mcsDist(1,2)/2);
    0225 minus=V.mcsDist(1,:)-(V.mcsDist(1,2)/2);
    0226
    0227 for j=1:size(plus,2);%columns
    0228
    0229     % spatial average ensembles on either side of the nodes
    0230     A(zi).Comp.mcsBack(:,j)=nanmean(A(zi).Clean.bsvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0231     % count the non-nan values that were averaged
    0232     A(zi).Comp.mcsBackContrib(:,j)=sum(~isnan(A(zi).Clean.bsvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0233
    0234     A(zi).Comp.mcsDir(:,j)=nanmean(A(zi).Clean.vDirvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0235     A(zi).Comp.mcsDirContrib(:,j)=sum(~isnan(A(zi).Clean.vDirvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0236
    0237     A(zi).Comp.mcsMag(:,j)=nanmean(A(zi).Clean.vMagvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0238     A(zi).Comp.mcsMagContrib(:,j)=sum(~isnan(A(zi).Clean.vMagvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0239
    0240     A(zi).Comp.mcsEast(:,j)=nanmean(A(zi).Clean.vEastvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0241     A(zi).Comp.mcsEastContrib(:,j)=sum(~isnan(A(zi).Clean.vEastvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0242
    0243     A(zi).Comp.mcsNorth(:,j)=nanmean(A(zi).Clean.vNorthvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0244     A(zi).Comp.mcsNorthContrib(:,j)=sum(~isnan(A(zi).Clean.vNorthvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0245
    0246     A(zi).Comp.mcsVert(:,j)=nanmean(A(zi).Clean.vVertvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0247     A(zi).Comp.mcsVertContrib(:,j)=sum(~isnan(A(zi).Clean.vVertvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0248
    0249     A(zi).Comp.mcsBed(:,j)=nanmean(A(zi).Clean.depthvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j))),2);
    0250     A(zi).Comp.mcsBedContrib(:,j)=sum(~isnan(A(zi).Clean.depthvecmap(:,(A(zi).Comp.itDist(1,:)<plus(1,j)&A(zi).Comp.itDist(1,:)>minus(1,j)))),2);
    0251
    0252 end
    0253
    0254 A(zi).Comp.mcsBack(A(zi).Comp.mcsBack>=255) = NaN;
    0255
    0256 end
    0257
    0258 %% Interpolate individual transects onto uniform mean c-s grid
    0259 % Fill in uniform grid based on individual transects mapped onto the mean
    0260 % cross-section by interpolating between adjacent points
    0261
    0262 %ZI = interp2(X,Y,Z,XI,YI)
    0263 for zi = 1 : z
    0264
    0265     A(zi).Comp.mcsBackI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0266         A(zi).Clean.bs(:,A(zi).Comp.vecmap),V.mcsDist, V.mcsDepth);
    0267     A(zi).Comp.mcsBackI(A(zi).Comp.mcsBackI>=255) = NaN;
    0268
    0269     A(zi).Comp.mcsDirI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0270         A(zi).Clean.vDir(:,A(zi).Comp.vecmap), V.mcsDist, V.mcsDepth);
    0271     A(zi).Comp.mcsMagI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0272         A(zi).Clean.vMag(:,A(zi).Comp.vecmap), V.mcsDist, V.mcsDepth);
    0273     A(zi).Comp.mcsEastI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0274         A(zi).Clean.vEast(:,A(zi).Comp.vecmap), V.mcsDist, V.mcsDepth);
    0275     A(zi).Comp.mcsNorthI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0276         A(zi).Clean.vNorth(:,A(zi).Comp.vecmap), V.mcsDist, V.mcsDepth);
    0277     A(zi).Comp.mcsVertI = interp2(A(zi).Comp.itDist, A(zi).Comp.itDepth, ...
    0278         A(zi).Clean.vVert(:,A(zi).Comp.vecmap), V.mcsDist, V.mcsDepth);
    0279
    0280     A(zi).Comp.mcsBedI  = interp1(A(zi).Comp.itDist(1,:),...
    0281         nanmean(A(zi).Nav.depth(A(zi).Comp.vecmap,:),2),V.mcsDist(1,:));
    0282
    0283 end
    0284
    0285 clear zi