function mkrawcdf(settings, outFileRoot)
% mkrawcdf - Read binary sonar logger file, write raw netcdf version
% mkrawcdf(metaFile, outFileRoot)
%
%(including the azimuth drive).  Outputs a netcdf file of the raw data
%
% Usage:
%    Input: settings are from read_globalatts and additions from the
%           sonar_meta.m file   
%              The program currently assumes azimuth,fan and pencil are 
%              processed as separate steps, so each needs it's own meta
%              file
%           outFileRoot - name of output netCDF file,
%              surrounded by single quotes WITHOUT the file
%              extension .cdf.  The files will be named:
%              'outFileRoot'_raw.cdf and 'outFileRoot _proc.cdf)
%
% Written by E.Montgomery, revisions by C. Sherwood
% Based on sonar processing software by C. Sherwood and C. Sullivan
% USGS Woods Hole Science Center
% emontgomery@usgs.gov
%
% Dependencies:
%   USGS NetCDF Toolbox (C. Denham)
%   readfan.m, readpencil.m, readrange.m  (M. Martini)
%   definerawsonarncfile.m  (E. Montgomery)
close all
more off

% get the current SVN version- the value is automatically obtained in svn
% is the file's svn.keywords is set to "Revision"
rev_info= 'SVN $Revision$';

% Check that the metadata contains required fields.  If a required field
% is missing, ask the user for it.
reqFields = {'FirstSonarDay','LastSonarDay','RootDataDir',...
   'SonartoAnimate'};
for f = 1:length(reqFields)
   if ~isfield(settings,reqFields{f})
      disp(['The field ''',reqFields{f},''' is not specified in ',metaFile,'.txt'])
      missingFields(f) = 1;
   else
      missingFields(f) = 0;
   end
end
if any(missingFields)
   disp('Required fields missing from the metadata-');
   disp (' check the metadata input file and try again')
   return;
   % Default settings for Eurostrat first deployment 2002-2003
end
clear reqFields missingFields

% Get a listing of sonar data directories
%expects directories of files by day under root data dir- outher directory
%structures aren't accommdated
SubDataDir = dir(settings.RootDataDir);

% Individual movie frames saved to a subdirectory
% of the current directory named 'frame'
outdir = 'frame';
[success,message,messageid] = mkdir(outdir);

clear FanTime PencilTime AzmTime
Fanidx = 1;
Pencilidx = 1;
Azmidx=1;
Ax1=1;
save settings settings;
% Go through each sonar data directory and parse. The first two
%  directories are always . and ..
for d=3:length(SubDataDir),
   %for d=3,
   disp(SubDataDir(d).name)

   % Get the year in which the data was collected
   [yyyy] = getSonarYear(SubDataDir(d).name, settings.FirstSonarDay,...
      settings.LastSonarDay);
   yyyy = str2num(yyyy);
   mm = str2num(SubDataDir(d).name(2:3));
   dd = str2num(SubDataDir(d).name(4:5));
   dirDay = julian(yyyy,mm,dd,0);
   fileDay = gregorian(dirDay - 1); %files are 1 day behind directory day
   settings.datayear = fileDay(1);
   fcnt=0;

   % Act only on directories
   if exist(fullfile(settings.RootDataDir,SubDataDir(d).name)) == 7,
      % Get a listing of data files
      DataFiles = dir(fullfile(settings.RootDataDir,SubDataDir(d).name));
      for f = 3:length(DataFiles),
         sonarDataFile = fullfile(settings.RootDataDir,...
            SubDataDir(d).name,...
            DataFiles(f).name);
         %disp(sonarDataFile)

         % get rid of all plots every 10th
         if mod(f,10)== 0; close all; end

         % The file name convention for sonar data files is
         % S1111020.F50 = S1MMDDHH.?MM where ? = F for fan
         % or P for pencil
         a = sscanf(DataFiles(f).name,'S1%c%c%c%c%c%c.%c%c%c');
         mon = str2num(a(1:2));
         day = str2num(a(3:4));
         hr = str2num(a(5:6));
         mn = str2num(a(8:9));
         datenum_sonar = datenum(settings.datayear,mon,day,hr,mn,0); 
        if (Ax1)     %only the first time through
             old_sonardate=datenum_sonar;
             Ax1=0;      
        end

         % The file name convention for individual
         % sonar frames is yyyymmddTHHMMSS.png
         outf = ['B',datestr(datenum_sonar,30)];
         outfile = fullfile(outdir,outf);
         if strmatch('azm',settings.SonartoAnimate)
            outFileRoot = ['az',datestr(datenum_sonar,29)];
         end

         % Fan files less than 5 kb are duds and
         % should not be used in the animation
         [WhichSonar, validFile] = checkSonarFile(sonarDataFile);

         % For valid (>5 kb) fan files, display the sonar image
         % and write the image data to NetCDF
         if WhichSonar == 'F' && validFile == 1 && strcmpi(settings.SonartoAnimate,'fan'),

            % first, display the sonar image
            FanTime(Fanidx) = julian([settings.datayear mon day hr mn 0]);
            FanDatenum=datenum([settings.datayear mon day hr mn 0]);

            % use new generic readfan
            fname=fullfile(settings.RootDataDir,SubDataDir(d).name,DataFiles(f).name)
            %  readfan and plotfan are replacing showfan07
            %  use readfan and plotfan OR showfan07 for the moment
            [FanData, FanHeader, FanSwitches] = readfan(fname);
            [npts,nscans]=size(FanData.imagedata);
            %
            % create the netcdf file
            if Fanidx==1
               % initiaize NetCDF dimensions
               ncdims.npoints=npts; ncdims.nscans=nscans;
               ofraw=[outFileRoot '_raw'];
               ncr = defineRawSonarNcFile(ofraw, settings, ncdims, 'r');
               ncr.start_time = datestr(gregorian(FanTime(1)));
               ncr=addHeaderMeta(ncr,FanHeader, FanSwitches);
               ncr.NOTE = ['vector variables combined from scan headers; ',...
                  'names match Imagenex documentation of headers']
               % also NOTE:
               % the StepSize in the FanSwitches (from which the attributes come)
               % is always WRONG, so compute it from the HeadAngles
                  ncr.StepSize(:)=median(diff(FanHeader.HeadAngle(100:150)));
               % update the history
               hist = ncr.history(:);
               hist_new = ['Sonar data written to NetCDF by ', mfilename, ', ' rev_info,...
                  ', using Matlab ', version, '; ',hist];
               ncr.history = hist_new;
               ncr = write2rawNcFile(ncr, Fanidx, FanTime(Fanidx), FanData, FanHeader);
            else
               ncr = write2rawNcFile(ncr, Fanidx, FanTime(Fanidx), FanData, FanHeader);
            end

            Fanidx=Fanidx+1;
            % save fan_data.mat FanData FanTime

            clear imagedata Xplot Yplot Zi
         elseif WhichSonar == 'P' && validFile == 1 && strcmpi(settings.SonartoAnimate,'pen'),
            % first, display the sonar image
            PencilTime(Pencilidx) = julian([settings.datayear mon day hr mn 0]);
            settings.PencilTime = PencilTime(Pencilidx);
            fname=fullfile(settings.RootDataDir,SubDataDir(d).name,DataFiles(f).name);

            % replaced showpencil07 with (readpencil and plotPencil07
            [PenData, PenHeader, PenSwitches] = readpencil(fname);

            % create the netcdf file
            if Pencilidx==1
               %  initiaize NetCDF attributes
               [npts,nscans]=size(PenData.imagedata);   % doesn't change size
               ncdims.npoints=npts; ncdims.nscans=nscans;
               ofraw=[outFileRoot '_raw'];
               % first do the raw
               ncr = defineRawSonarNcFile(ofraw, settings, ncdims,'r');
               ncr.start_time = datestr(gregorian(PencilTime(1)));
                  ncr=addHeaderMeta(ncr,PenHeader, PenSwitches);
               hist = ncr.history(:);
               hist_new = ['Sonar data written to NetCDF by ', mfilename, ', ' rev_info,...
                  ', using Matlab ', version, '; ',hist];
               ncr.history = hist_new;
               %  don't know why, but the flipud is needed here
               ncr.NOTE = ['Pencil data']
               % note StepSize in PenSwitches is CORRECT, so no
               % computation is necessary here.
               ncr = write2rawNcFile(ncr, Pencilidx, PencilTime(Pencilidx), PenData, PenHeader);
            else
               ncr = write2rawNcFile(ncr, Pencilidx, PencilTime(Pencilidx),PenData, PenHeader);
            end
            Pencilidx = Pencilidx+1;
            if Pencilidx>=293    % for MVCO-07, this is after the last good data
               disp ('Pencil idx ge 293')
            end

            clear outfile imfile
            % now add a section to do the Azimut drive data
         elseif WhichSonar == 'A' && validFile == 1 && strcmpi(settings.SonartoAnimate,'azm')
           %these files bomb if get too big, so forcing into 1 day per file.
           % check to see if day has changed BEFORE reading new record
            if floor(datenum_sonar)==floor(old_sonardate);
               clear outfile AzmData AzmHeader AzmSwitches
            else
               %close existing .ncfile
                  %end for AzmTime should be 4, eventhough Azmidx is 5
               ncr.stop_time = datestr(gregorian(AzmTime(end)));
               if isempty(ncr.DELTA_T(:))
                   ncr.DELTA_T = [num2str(gmedian(diff(AzmTime))*60*60*24),' sec'];
               end
               ncr.azm_step= [num2str(gmean(diff(AzmHeader.AAngle(:,:))))];
               close(ncr);
               % make it so there'll be a new file with date as part of it
               outFileRoot = ['az',datestr(datenum_sonar,29)];
               Azmidx=1;
               % AzmTime should only the times for this day, so clear each close
               clear AzmTime
            end
            % Azm time is derived from the filename
            AzmTime(Azmidx) = julian([settings.datayear mon day hr mn 0]);
            settings.StartRecordTime = datestr(gregorian(AzmTime(Azmidx)));
            settings.AzmTime = AzmTime(Azmidx);
            fname=fullfile(settings.RootDataDir,SubDataDir(d).name,DataFiles(f).name);

            % use readrangeall to get one set of azimuth sonar sweeps
            mkplt='n';
            [AzmData, AzmHeader, AzmSwitches] = readrangeall(fname, mkplt);

            % create the netcdf file
            if Azmidx==1
               %  initiaize NetCDF attributes
               [nrts,npts,nscans]=size(AzmData.imagedata);   % doesn't change size
               ncdims.rots=nrts; ncdims.npoints=npts; ncdims.nscans=nscans;
               ofraw=[outFileRoot '_raw'];
               % initialize the nc file and write some metadata
               ncr = defineRawSonarNcFile(ofraw, settings, ncdims,'r');
               ncr.start_time = datestr(gregorian(AzmTime(Azmidx)));
               ncr=addHeaderMeta(ncr,AzmHeader, AzmSwitches);
		       hist = ncr.history(:);
               hist_new = ['Sonar data written to NetCDF by ', mfilename, ', ' rev_info,...
                  ', using Matlab ', version, '; ',hist];
               ncr.history = hist_new;
               ncr.NOTE = ['Azimuth data'];
               % note StepSize in AzmSwitches is CORRECT, so no
               % computation is necessary here.
               ncr = write2rawNcFile(ncr, Azmidx, AzmTime(Azmidx), AzmData, AzmHeader);
            else
               ncr = write2rawNcFile(ncr, Azmidx, AzmTime(Azmidx),AzmData, AzmHeader);
            end
            Azmidx=Azmidx+1;
            old_sonardate = datenum_sonar;
         end
      end
   end
end

% this is where the data is saved
if strcmpi(settings.SonartoAnimate,'fan'),
   %   save fan_data.mat FanData FanTime FanHeader
   ncr.stop_time = datestr(gregorian(FanTime(end)));
   if isempty(ncr.DELTA_T(:))
       ncr.DELTA_T = [num2str(gmean(diff(FanTime*1000))*60),' sec'];
   end
   close(ncr);
   ncclose;
elseif strcmpi(settings.SonartoAnimate,'pen'),
   %   save pen_data.mat PencilTime PenData PenHeader
   ncr.stop_time = datestr(gregorian(PencilTime(end)));
   if isempty(ncr.DELTA_T(:))
     ncr.DELTA_T = [num2str(gmean(diff(PencilTime*1000))*60),' sec'];
   end
   close(ncr);
   ncclose;
% azimuth file will probably be closed by the time it gets here   
elseif strcmpi(settings.SonartoAnimate,'azm'),
    % see if the file's still open - sz(1) will be >0, if the nc file is open
    sz=size(ncr);
    if sz(1)>1
        %   save azm_data.mat AzmTime AzmData AzmHeader
        ncr.stop_time = datestr(gregorian(AzmTime(end)));
        if isempty(ncr.DELTA_T(:))
            ncr.DELTA_T = [num2str(gmean(diff(AzmTime*1000))*60),' sec'];
        end
        close(ncr);
        ncclose;
    end
end
% ---------------- Subfunction: readSonarMeta.m ------------------------- %
function userMeta = readSonarMeta(metaFile)
[atts, defs] = textread(metaFile,'%s %63c','commentstyle','shell');
defs = cellstr(defs);
for i = 1:length(atts)
   theAtt = atts{i}(:)';
   theDef = defs{i}(:)';
   % deblank removes trailing whitespace
   theAtt = deblank(theAtt);
   theDef = deblank(theDef);
   % check for and replace spaces in
   % the attributes with underscores
   f1 = find(isspace(theAtt));
   f2 = strfind(theAtt,'-');
   f = union(f1,f2);
   if ~isempty(f)
      theAtt(f) = '_';
   end
   % attribute definitions read in as characters; convert to
   % numbers where appropriate
   theDefNum = str2double(theDef);
   if ~isnan(theDefNum) %if the conversion worked (field was numeric)
       if ~strcmp(theAtt,'MOORING');  % mooring should remain a string
         theDef = theDefNum;
       end
   end
   eval(['userMeta.',theAtt,'= theDef;'])
end

% ---------------- Subfunction:  getSonarYear.m ------------------------- %
function [yyyy] = getSonarYear(DirName, FirstSonarDay, LastSonarDay)

first = datevec(FirstSonarDay);
first_yyyy = first(1);
first_mm = first(2);
first_dd = first(3);

last = datevec(LastSonarDay);
last_yyyy = last(1);
last_mm = last(2);
last_dd = last(3);

mm = str2num(DirName(2:3));
dd = str2num(DirName(4:5));

if (mm >= first_mm)
   yyyy = num2str(first_yyyy);
else
   yyyy = num2str(last_yyyy);
end

% ---------------- Subfunction:  checkSonarFile.m ----------------- %
function [WhichSonar, validFile] = checkSonarFile(sonarFileName)

F = dir(sonarFileName); % FAN files only!
filesize = F.bytes;
if strcmp(F.name(10),'F') && filesize <= 5000
   disp(['### The Fan file ',sonarFileName, ' is a dud file and will not be used in animations'])
   WhichSonar = 'F';
   validFile = 0;
elseif strcmp(F.name(10),'F') && filesize > 5000
   WhichSonar = 'F';
   validFile = 1;
elseif strcmp(F.name(10),'P') && filesize > 100
   WhichSonar = 'P';
   validFile = 1;
elseif strcmp(F.name(10),'R') && filesize > 100
   WhichSonar = 'A';
   validFile = 1;
else
   WhichSonar = 'N'; %  N== NG need to pass back something here
   validFile = 0;
end
% ---------------- Subfunction:  addHeaderMeta.m -------------------- %
function nc =addHeaderMeta(nc, header, switches)
% add everything from the switches to the global attributs
nc.HeadID=switches.HeadID;
nc.Range=switches.Range;
nc.RevHold=ncchar(switches.RevHold);
nc.MasterSlave=ncchar(switches.MasterSlave);
nc.StartGain=switches.StartGain;
nc.LOGF=switches.LOGF;
nc.Absorption=switches.Absorption;
nc.TrainAngle=switches.TrainAngle;
nc.SectorWidth=switches.SectorWidth;
nc.StepSize=switches.StepSize;
nc.PulseLength=switches.PulseLength;
nc.DataPoints=switches.DataPoints;
nc.DataBits=switches.DataBits;
nc.UpBaud=switches.UpBaud;
nc.Profile=ncchar(switches.Profile);
nc.Calibrate=ncchar(switches.Calibrate);
nc.SwitchDelay=switches.SwitchDelay;
% these are fan only
if isfield(header,'FanTime')
    nc.RangeOffset=switches.RangeOffset;
    nc.StepDirection=switches.StepDirection;
    nc.MoveRelative=switches.MoveRelative;
    % frequency is not returned in the fan switches
    nc.Frequency=ncchar(['675 kHz']);
end
% azimuth doesn't have these, pen and fan do.
if isempty(strfind(header.FileName,'.R'))
    cc=char(header.ReturnDataHeaderType{1});
    nc.HeadType=ncchar(cc);
    nc.NDataBytes=header.NDataBytes(1);
    nc.NReturnBytes=header.NReturnBytes(1);
end
% only pencil and azimuth have frequency setting returned
if isempty(strfind(header.FileName,'.F'))
    nc.Frequency=ncchar([num2str(switches.Frequency) ' kHz']);
end
% ---------------- Subfunction:  write2rawNcFile.m -------------------- %
function nc = write2rawNcFile(nc, idx, timeWrd, data, header)
% treat fan and pencil the same.

% timeWrd is in decimal days, so have to convert .604 days to msecs
nc{'time'}(idx) = floor(timeWrd);
nc{'time2'}(idx) = mod(timeWrd,floor(timeWrd))*1000*3600*24;
nn=size(data.imagedata);
if length(nn)==2
   [npoints,nscans]=size(data.imagedata);
else
   [nrot,npoints,nscans]=size(data.imagedata);
end
if idx ==1
   nc{'points'}(1:npoints)=[1:1:npoints];
   nc{'scan'}(1:nscans)=[1:1:nscans];
   if length(nn) == 3
      nc{'rots'}(1:nrot)=[1:1:nrot];
   end
end
if length(nn) == 2
   if idx==1
      nc{'headangle'}(1:nscans) = header.HeadAngle;
      nc{'nDataBytes'}(1:nscans) = header.NDataBytes;
      nc{'head_pos'}(idx,1:nscans) = header.HeadPosition;
      nc{'profile_range'}(idx,1:nscans) = header.ProfileRange;
      nc{'raw_image'}(idx,1:npoints,1:nscans) = data.imagedata;
   else
      nc{'head_pos'}(idx,1:nscans) = header.HeadPosition;
      nc{'profile_range'}(idx,1:nscans) = header.ProfileRange;
      nc{'raw_image'}(idx,1:npoints,1:nscans) = data.imagedata;
   end
elseif length(nn) == 3 %azimuth
   nc{'headangle'}(idx,1:nrot,1:nscans-1) = header.HeadAngle;
   nc{'profile_range'}(idx,1:nrot,1:nscans-1) = header.ProfileRange;
   nc{'azangle'}(idx,1:nrot) = header.AAngle;
   nc{'raw_image'}(idx,1:nrot,1:npoints,1:nscans) = data.imagedata;
end

clear alltime timeWrd