% metaADVExample - an example script for processing a single ADV data file.  
% The steps here follow the steps outlined in the manual
%
% you should be running this script from the directory in which the data
% files reside

% Written by Marinna Martini for the U.S. Geological Survey 
% Coastal & Marine Program Woods Hole Field Center, Woods Hole, MA
% http://woodshole.er.usgs.gov/ Please report bugs to mmartini@usgs.gov

% ------------- Step I what are the file names?
% processing file for 7772001.adr
adrFile = '7772001.adr';
outFileRoot = '7772adv'; % the name for the output files

% ------------- Step II what is the metadata?
%
% here is the metadata structure for this deployment
%  - if you use an empty char string denoted by '' for an external sensor serial
%    number, then that data will not be included
%  - calibration info will be assigned to variable attributes dynamically
%    based on the field names of items under the struct field cals.
% Do not nest more structs under cals, or adr2cdf will not record the data
%
% no two data loggers can have the same four digit mooring number
metadata.mooring_number = '7772'; %four digit NNNL, NNN mooring + L logger
metadata.deployment_start = '09-feb-05 17:00:00'; % the in water time
metadata.deployment_end = '18-may-2005 13:02 '; % the out of water time
metadata.lon = -70.7809; % decimal degrees, West = negative
metadata.lat = 42.3786; % decimal degrees, South = negative
metadata.declination = -(15+35/60); % degrees west is negative
metadata.water_depth = 104.9*0.3048; % m (DO NOT OMIT!)
metadata.system.serial = 'G345'; % serial number of the Hydra electronics
metadata.system.synch_to = 'G357'; % G### | nothing
metadata.adv.serial = 'B227';
metadata.adv.height = 1.203; % height of emitter above bottom, m
                            % used to reference anything coming from the ADV probe
% distance from emitter to sample volume: 18 cm for ADVO, 10 cm for ADVF 10 MHz, 5 cm for ADVF 16 MHz                           
metadata.adv.sample_volume_offset = .18; % m
% for attached sensors, use '' for the serial number if there is no sensor, 
% and variable will not be written
metadata.extpress.serial = 'P68022'; % use the druck file name without the extension
metadata.extpress.height = 2.14; % m, same as adv if in probe
% pressure sensor calibration data will be read from the .drk file
% now add the external sensors, in this example, the OBS was on channel 1
% and the trans on channel 2
% anything under the metadat.ext#.cals structure gets written as an attibute
% in the netCDF file
metadata.ext1.serial = 'TRANS289'; % OBS### | TRANS### | XXX####
metadata.ext1.height = 1.81; % m
metadata.ext1.cals.manufacturer = 'SeaTech'; % manufacturer
% the following are required for trans data
metadata.ext1.cals.path_length = 0.25; % focal length of the trans, 0.5 | 0.25 | 1 m, required!
metadata.ext1.cals.pre = [4.465 0]; % for trans, air reading: [unblocked blocked]
metadata.ext1.cals.post = [4.492 0]; % for trans, air reading: [unblocked blocked]
% anything under the metadat.ext#.cals structure gets written as an attibute
% in the netCDF file
metadata.ext2.serial = 'OBS2183'; % OBS### | TRANS### | XXX####
metadata.ext2.height = 0.98; % m
metadata.ext2.cals.gain = .33; % for OBS gain setting, g/l per volt
% if coef is present, dohydraoptic will apply the coefficients to OBS data
% use [0 1 0], [] or omit if there is no calibration for concentration
metadata.ext2.cals.coef = []; % coefficients to a calibration eqation
%metadata.ext2.cals.units = 'kg/m-3'; % or g/l
%metadata.ext2.cals.equation = 'Conc [kg/m-3] = coef(1)+coef(2)*V+coef(3)*V^2';
% the following are required for trans data
metadata.CTD.serial = ''; % 
metadata.CTD.height = 0; % m
% SBE37's record data in real units and do not need applied cals
metadata.metafile = mfilename('fullpath');
metadata.metafile_date = datestr(now);
metadata.metafile_version = '0.0';
metadata.metafile_author = 'MM/EM';
% if the following fields are missing, UNKOWN will be used
metadata.origin = 'USGS/WHFC'; % with collaborator's data, could be USC, etc.
metadata.experiment = 'Boston Longterm'; 
metadata.project = 'WHFC'; % might also use OFA funding agency, such as MWRA, EPA, WCMG
metadata.whichscheme = 1; % which sampling scheme to read out
% always 1 for single sampling scheme data, not more than 3 

% it is useful to omit bursts at the beginning and the end that were not recorded
% when the instrument was in the water
% to process all data use : metadata.goodbursts = [1 Inf];
metadata.goodbursts = [207 2839];

% these programs take a long time to run, so keep a diary
diary(sprintf('run%s',datestr(now,30))) 

% change the 0's to 1's sequentially for the steps below to work through
% the processing.  Most of the programs take several hours to run on 30K bursts, so
% testing the threshold's used on a small subset of the data is encouraged
% before processing everything.

% ---------- Step III convert from vendor's raw binary to netCDF & get quality info
if 0,
    % now run the file to convert from raw to netcdf
    diagnostics.pcountout = 0;      % both should be 0 !!!
    diagnostics.fixpfreq = 0;
    adr2cdf(adrFile,outFileRoot,metadata.goodbursts,metadata.goodbursts,diagnostics,metadata);
end

% ---------- Step IV automated cleaning
if 0,
    % make a copy of the files, fixbadadv works directly on them
    % use the already canged *s1.cdf to be changed further
    disp(sprintf('copying %s to %s',[outFileRoot 'b.cdf'], [outFileRoot 'b1.cdf']))
    copyfile([outFileRoot 'b.cdf'], [outFileRoot 'b1.cdf'])
    disp(sprintf('copying %s to %s',[outFileRoot 's.cdf'], [outFileRoot 's1.cdf']))
    copyfile([outFileRoot 's.cdf'], [outFileRoot 's1.cdf'])
    disp(sprintf('copying %s to %s',[outFileRoot 'q.cdf'], [outFileRoot 'q1.cdf']))
    copyfile([outFileRoot 'q.cdf'], [outFileRoot 'q1.cdf'])
end

if 0,
    % function cleanhydra(cdfbFile, cdfsFile, cdfqFile, operation, ...
    %       'variables', {vnames}, 'settings', settings.n, ...
    %       'burstrange', [1:Inf]) 
    % the process below is reasonably robust over multiple data sets
    % to run thumbfinger with std of 10, use this:
      seta.nsd=10;
      seta.nsd_pr=20;   % set to a large number to do nothing to good pressures
       cleanhydra([outFileRoot 'b1.cdf'], [outFileRoot 's1.cdf'], ...
        [outFileRoot 'q1.cdf'], 'thumbfinger','settings',seta);
    % make a copy of the files, fixbadadv works directly on them
    % use the already canged *s1.cdf to be changed further
    disp(sprintf('copying %s to %s',[outFileRoot 'b1.cdf'], [outFileRoot 'b2.cdf']))
    copyfile([outFileRoot 'b1.cdf'], [outFileRoot 'b2.cdf'])
    disp(sprintf('copying %s to %s',[outFileRoot 's1.cdf'], [outFileRoot 's2.cdf']))
    copyfile([outFileRoot 's1.cdf'], [outFileRoot 's2.cdf'])
    disp(sprintf('copying %s to %s',[outFileRoot 'q1.cdf'], [outFileRoot 'q2.cdf']))
    copyfile([outFileRoot 'q1.cdf'], [outFileRoot 'q2.cdf'])

    % to run just deglitch1vector with nsd=2.8, use:
      seta.nsd=2.8;
       cleanhydra([outFileRoot 'b2.cdf'], [outFileRoot 's2.cdf'], ...
        [outFileRoot 'q2.cdf'], 'deglitch1vector','settings', seta);
 
    % to modify the default settings, and just do bursts 1361-1365, use this syntax:
    %  seta.nsd=3; seta.ndt=10;
    %   cleanhydra([outFileRoot 'b1.cdf'], [outFileRoot 's1.cdf'], [outFileRoot 'q11.cdf'], ...
    %    'deglitch1vector','settings',seta,'burstrange', [1361:1365]);

    % use of autoclean (does thumbfinger, corrcheck, then deglitch all using nsd=2.8)
    %  is possible, but not recommended because the results are awful
    %   cleanhydra([outFileRoot 'b1.cdf'], [outFileRoot 's1.cdf'], ...
    %    [outFileRoot 'q1.cdf'], 'autoclean');
    
end

% ------------ Step V interactive flagging
%
%  This instrument got +/- 100 vels in storms, which have very high std in
%  those periods, but the burst data looks fine-  scutoff of 40, 55 and 9
%  were used to allow all bursts.
%  ** for other datasets lower thresholds are preferable **
if 0,
    % settings.currentBeam is 1 by default
    seta.currentBeam=1;   seta.new_scutoff= [40 40 40 40 40]; % match values from beam 1
      flagbadadv('QualityFile', [outFileRoot 'q2.cdf'], 'settings',seta)
    seta.currentBeam=2;   seta.new_scutoff= [55 55 55 55 55]; % match values from beam 1
      flagbadadv('QualityFile', [outFileRoot 'q2.cdf'], 'settings',seta)
    seta.currentBeam=3;   seta.new_scutoff= [9 9 9 9 9]; % match values from beam 1 & 2
      flagbadadv('QualityFile', [outFileRoot 'q2.cdf'], 'settings',seta)
end

if 0,
    % fixbadadv(dirtybFile, dirtysFile, Qdata)
    fixbadadvvel([outFileRoot 'b2.cdf'],[outFileRoot 's2.cdf'],[outFileRoot 'q2.cdf']);
end

% ---------- Step VI convert to scientific units and apply rotations
if 0,
    % another time consuming step, so keep a diary
    diary(sprintf('run%s',datestr(now,30)))
    % unless otherwise set...
    % switches.ecorr_showplots = 0;
    % deglitch has been removed... use cleanhydra
    % despike has been removed... use cleanhydra
    % orientation overrides left intact in adv2nc
    % Vector length must match the number of bursts this function is
    %  asked to process
    % switches.override_heading = [];
    % switches.override_pitch = [];
    % switches.override_roll = [];
    % pressure mask left intact in adv2nc
    % switches.override_pmask = []; % pressure spike mask
    % example to omit the first 5 samples: 
    % the ones vector must be sized to match the number of samples in the burst
    %    switches.override_pmask = [zeros(5,1); ones(480-5,1)];
    % to remove the minimum voltage from obs prior to do_hydraoptic, set
    % this to 1
    % switches.remove_minVolt=0;
    % switches = []; % invoke all defaults

    Adv = adv2nc([outFileRoot 'b2.cdf'], [outFileRoot 's2.cdf'], outFileRoot, switches, [207 28339]);
% this run took 415 minutes, so don't expect fast completion!
end

diary off