Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Bug binmapping #777

Merged
merged 2 commits into from
Apr 13, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
60 changes: 60 additions & 0 deletions Preprocessing/+NortekADCP/adcp_adjusted_distances_nortek3beam.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,60 @@
function nonMappedHeightAboveSensorBeam = adcp_adjusted_distances_nortek3beam(roll, pitch, distAlongBeams, beamAngle, number_of_beams)
% function nonMappedHeightAboveSensorBeam = adcp_adjusted_distances_nortek3beam(roll, pitch, distAlongBeams, beamAngle, number_of_beams)
%
% Adjusts height above sensor for each beam in adcpBinMapping.m
% Can read 3 beam Nortek instruments.
% Uses the beam angle and the roll/pitch of the instrument.
% No adjustement when tilt is zero,
% nonMappedHeightAboveSensorBeam will be the same than distAlongBeams.
%
% Inputs:
%
% roll - 1D single precision array (Ntime), the internal roll angles in degrees
% pitch - 1D single precision array (Ntime), the internal pitch angles in degrees
% distAlongBeams - 1D array (Nbins), the internal distance along beam for each bin
% beamAngle - scalar, beam angle from the instrument
% number_of_beams - interger, number of beams on the ADCP
%
% Outputs:
%
% nonMappedHeightAboveSensorBeam - 3D single precision array (Ntime, Nbins, Nbeams),
% tilt adjusted height above sensor for
% each of the 3 beams
%
% Example:
%
% roll_ini = [-180:10:180]';
% pitch_ini = [-180:10:180]';
% distAlongBeams = [20:20:600];
% beamAngle = 0.3;
% number_of_beams = 3;
%
% x = NortekADCP.adcp_adjusted_distances_nortek3beam(roll_ini,pitch_ini,distAlongBeams,beamAngle,number_of_beams);
%
% idx0 = find(roll_ini == 0);
% for i=1:number_of_beams; assert(isequal(x(idx0,:,i), distAlongBeams)); assert(~isequal(x(~idx0,:,i), distAlongBeams)); end
%
%

narginchk(5,5)

number_of_beams = 3;
nBins = length(distAlongBeams);

% set to single precision
nonMappedHeightAboveSensorBeam = nan(length(pitch), length(distAlongBeams), number_of_beams, 'single');

nonMappedHeightAboveSensorBeam(:,:,1) = (cos(beamAngle - pitch) / cos(beamAngle)) * distAlongBeams;
nonMappedHeightAboveSensorBeam(:,:,1) = repmat(cos(roll), 1, nBins) .* nonMappedHeightAboveSensorBeam(:,:,1);

beamAngleX = atan(tan(beamAngle) * cos(60 * pi / 180)); % beams 2 and 3 angle projected on the X axis
beamAngleY = atan(tan(beamAngle) * cos(30 * pi / 180)); % beams 2 and 3 angle projected on the Y axis

nonMappedHeightAboveSensorBeam(:,:,2) = (cos(beamAngleX + pitch) / cos(beamAngleX)) * distAlongBeams;
nonMappedHeightAboveSensorBeam(:,:,2) = repmat(cos(beamAngleY + roll) / cos(beamAngleY), 1, nBins) .* nonMappedHeightAboveSensorBeam(:,:,2);

nonMappedHeightAboveSensorBeam(:,:,3) = (cos(beamAngleX + pitch) / cos(beamAngleX)) * distAlongBeams;
nonMappedHeightAboveSensorBeam(:,:,3) = repmat(cos(beamAngleY - roll) / cos(beamAngleY), 1, nBins) .* nonMappedHeightAboveSensorBeam(:,:,3);

end

75 changes: 75 additions & 0 deletions Preprocessing/+TeledyneADCP/adcp_adjusted_distances_rdi4beam.m
Original file line number Diff line number Diff line change
@@ -0,0 +1,75 @@
function nonMappedHeightAboveSensorBeam = adcp_adjusted_distances_rdi4beam(roll, pitch, pitchSign, distAlongBeams, beamAngle, number_of_beams)
% function nonMappedHeightAboveSensorBeam = adcp_adjusted_distances_rdi4beam(roll, pitch, pitchSign, distAlongBeams, beamAngle, number_of_beams)
%
% Adjusts height above sensor for each beam in adcpBinMapping.m
% Can read 4 beam RDI instruments.
% Uses the beam angle and the roll/pitch of the instrument.
% No adjustement when tilt is zero,
% nonMappedHeightAboveSensorBeam will be the same than distAlongBeams.
% Deals with up and down facing ADCPs
%
% Inputs:
%
% roll - 1D single precision array (Ntime), the internal roll angles in degrees
% pitch - 1D single precision array (Ntime), the internal pitch angles in degrees
% pitchSign - integer (-1 or 1), pitch sign based on the ADCP's orientation
% distAlongBeams - 1D array (Nbins), the internal distance along beam for each bin
% beamAngle - scalar, beam angle from the instrument
% number_of_beams - interger, number of beams on the ADCP
%
% Outputs:
%
% nonMappedHeightAboveSensorBeam - 3D single precision array (Ntime, Nbins, Nbeams),
% tilt adjusted height above sensor for
% each of the 4 beams
%
% Example:
%
% %testing for up and down facing ADCPs - no change for beam 1 and 2 (roll)
% %changes only for beam 3 and 4 (pitch)
%
% roll_ini = [-180:10:180]';
% pitch_ini = [-180:10:180]';
% pitchSign_up = 1;
% pitchSign_down = -1;
% distAlongBeams = [20:20:600];
% beamAngle = 0.3;
% number_of_beams = 4;
%
% x_up = TeledyneADCP.adcp_adjusted_distances_rdi4beam(roll_ini,pitch_ini,pitchSign_up,distAlongBeams,beamAngle,number_of_beams);
% x_down = TeledyneADCP.adcp_adjusted_distances_rdi4beam(roll_ini,pitch_ini,pitchSign_down,distAlongBeams,beamAngle,number_of_beams);
%
% assert(isequal(x_up(:,:,1), x_down(:,:,1)));
% assert(isequal(x_up(:,:,2), x_down(:,:,2)));
% assert(~isequal(x_up(:,:,3), x_down(:,:,3)));
% assert(~isequal(x_up(:,:,4), x_down(:,:,4)));
%
% idx0 = find(roll_ini == 0);
% for i=1:number_of_beams; assert(isequal(x_up(idx0,:,i), distAlongBeams)); assert(isequal(x_down(idx0,:,i), distAlongBeams)); end
% for i=1:number_of_beams; assert(~isequal(x_up(~idx0,:,i), distAlongBeams)); assert(~isequal(x_down(~idx0,:,i), distAlongBeams)); end
%
%

narginchk(6,6)

nmh=1;
if pitchSign < 0
% adjust nonMappedHeightAboveSensorBeam sign when downfacing
nmh = -1;
% adjust roll when downfacing
roll = roll-pi;
end

CP = cos(pitchSign * pitch);
CR = cos(roll);

% set to single precision
nonMappedHeightAboveSensorBeam = nan(length(pitch), length(distAlongBeams), number_of_beams, 'single');

nonMappedHeightAboveSensorBeam(:,:,1) = nmh * (CP .* (cos(beamAngle + roll) / cos(beamAngle) .* distAlongBeams));
nonMappedHeightAboveSensorBeam(:,:,2) = nmh * (CP .* (cos(beamAngle - roll) / cos(beamAngle) .* distAlongBeams));
nonMappedHeightAboveSensorBeam(:,:,3) = nmh * (CR .* (cos(beamAngle - pitchSign * pitch) / cos(beamAngle) .* distAlongBeams));
nonMappedHeightAboveSensorBeam(:,:,4) = nmh * (CR .* (cos(beamAngle + pitchSign * pitch) / cos(beamAngle) .* distAlongBeams));

end

50 changes: 18 additions & 32 deletions Preprocessing/adcpBinMappingPP.m
Original file line number Diff line number Diff line change
Expand Up @@ -69,9 +69,9 @@
% RDI 4 beams ADCPs:
% It is assumed that the beams are in a convex configuration such as beam 1
% and 2 (respectively 3 and 4) are aligned on the pitch (respectively roll)
% axis. When pitch is positive beam 3 is closer to the surface while beam 4
% gets further away. When roll is positive beam 2 is closer to the surface
% while beam 1 gets further away.
% axis. For an upward looking ADCP, when pitch is positive beam 3 is closer
% to the surface while beam 4 gets further away. When roll is positive beam
% 2 is closer to the surface while beam 1 gets further away.
%
% Nortek 3 beams ADCPs:
% It is assumed that the beams are in a convex configuration such as beam 1
Expand All @@ -81,47 +81,33 @@
% while beams 2 and 3 get further away. When roll is positive beam 3 is
% closer to the surface while beam 2 gets further away.
%

% Set pitch sign to deal with downward facing ADCPs
pitchSign = 1;

distAlongBeams = sample_data{k}.dimensions{distAlongBeamsIdx}.data';
if all(diff(distAlongBeams)<0)
%invert distAlongBeams so we have increasing values
% this is required for the interpolation function.
distAlongBeams = distAlongBeams*-1;
pitchSign = -1;
end
pitch = sample_data{k}.variables{pitchIdx}.data * pi / 180;
roll = sample_data{k}.variables{rollIdx}.data * pi / 180;
beamAngle = sample_data{k}.meta.beam_angle * pi / 180;

%TODO: the adjusted distances should be exported
% as variables to enable further diagnostic plots and debugging.
%TODO: the adjusted distances should be a Nx4 array or Nx3 array.
if isRDI% RDI 4 beams
if isRDI && absic4Idx ~= 0 % RDI 4 beams
number_of_beams = 4;
%TODO: block function.
CP = cos(pitch);
% H[TxB] = P[T] x (+-R[T] x B[B])
nonMappedHeightAboveSensorBeam1 = CP .* (cos(beamAngle + roll) / cos(beamAngle) .* distAlongBeams);
nonMappedHeightAboveSensorBeam2 = CP .* (cos(beamAngle - roll) / cos(beamAngle) .* distAlongBeams);

% H[TxB] = R[T] x (-+P[T] x B[B])
CR = cos(roll);
nonMappedHeightAboveSensorBeam3 = CR .* (cos(beamAngle - pitch) / cos(beamAngle) .* distAlongBeams);
nonMappedHeightAboveSensorBeam4 = CR .* (cos(beamAngle + pitch) / cos(beamAngle) .* distAlongBeams);
else
number_of_beams = 3;
nBins = length(distAlongBeams);
%TODO: block function, include tests.
% Nortek 3 beams
nonMappedHeightAboveSensorBeam1 = (cos(beamAngle - pitch) / cos(beamAngle)) * distAlongBeams;
nonMappedHeightAboveSensorBeam1 = repmat(cos(roll), 1, nBins) .* nonMappedHeightAboveSensorBeam1;

beamAngleX = atan(tan(beamAngle) * cos(60 * pi / 180)); % beams 2 and 3 angle projected on the X axis
beamAngleY = atan(tan(beamAngle) * cos(30 * pi / 180)); % beams 2 and 3 angle projected on the Y axis
nonMappedHeightAboveSensorBeam = TeledyneADCP.adcp_adjusted_distances_rdi4beam(roll, pitch, pitchSign, distAlongBeams, beamAngle, number_of_beams);

elseif isNortek && absic4Idx == 0 % Nortek 3 beams
number_of_beams = 3;

nonMappedHeightAboveSensorBeam2 = (cos(beamAngleX + pitch) / cos(beamAngleX)) * distAlongBeams;
nonMappedHeightAboveSensorBeam2 = repmat(cos(beamAngleY + roll) / cos(beamAngleY), 1, nBins) .* nonMappedHeightAboveSensorBeam2;
nonMappedHeightAboveSensorBeam = NortekADCP.adcp_adjusted_distances_nortek3beam(roll, pitch, distAlongBeams, beamAngle, number_of_beams);

nonMappedHeightAboveSensorBeam3 = (cos(beamAngleX + pitch) / cos(beamAngleX)) * distAlongBeams;
nonMappedHeightAboveSensorBeam3 = repmat(cos(beamAngleY - roll) / cos(beamAngleY), 1, nBins) .* nonMappedHeightAboveSensorBeam3;
end

nSamples = length(pitch);
Expand All @@ -136,13 +122,13 @@

switch n
case 1
dvar = nonMappedHeightAboveSensorBeam1;
dvar = nonMappedHeightAboveSensorBeam(:,:,n);
case 2
dvar = nonMappedHeightAboveSensorBeam2;
dvar = nonMappedHeightAboveSensorBeam(:,:,n);
case 3
dvar = nonMappedHeightAboveSensorBeam3;
dvar = nonMappedHeightAboveSensorBeam(:,:,n);
case 4
dvar = nonMappedHeightAboveSensorBeam4;
dvar = nonMappedHeightAboveSensorBeam(:,:,n);
otherwise
errormsg('Beam number %s not supported', num2str(n))
end
Expand Down