function gps2gbin(rawfilebase, varargin)
% gps2gbin(rawfilebase, varargin)
%                                
% read in a .raw file, read in the .raw.log file to get logger times,
% read in enough  .gps.bin files to get logger times and gga times
% correct the raw times to logger times 
% correct the logger times to gga times (if timebase == 'gga')
% - or - correct gga times to logger times (if timebase == 'log')
% write a flat binary file with 3xN entries
%   [filtered_loggercorrected_rawdday   bestdday sec_pc_minus_utc ]'
% output name is 'outfpath/rawfilebase.tim.gbin' for "gridded bin"
% 
% eg. gps2gbin('nbp2001_094_36014')
% varargin name, value pairs are:
% 
% binfpath [path to *.gps.bin]       %default 'binfiles/pcode'
% outfpath [path to output file]     %default 'gbinfiles/pcode'
% rawfpath [path to *.raw]           %default 'raw/adcp_pings'
% timebase ['gga' or 'log']          %default 'log'
%           use 'log' if the logger
%           is slaved to a timeserver
% yearbase                           %default (from rawfilebase)
% gap_seconds [NaN positions with    %default 10
%             gaps longer than this]
% scootch_sec [require time shifts to         % .5
%              be under this
% back_to_front [0 or 1]             %default is 0 (slower but mroe
%                                        reliable front_to_back)
% redo     [0 or 1] redo fies        %default 0
% instclass [nb]                     %default nb 
%           (not imlimented: os, bb, son)
% verbose [0 or 1]                   %default off

cfg.binfpath = 'binfiles/pcode';
cfg.outfpath = 'gbinfiles/time';
cfg.rawfpath = 'raw/adcp_pings';
cfg.timebase = 'log';
cfg.back_to_front = 0;
cfg.gap_seconds = 10;
cfg.scootch_sec = .5;
cfg.instclass = nb;
cfg.redo = 0;
cfg.verbose  = 0;
all_ = find('_' == rawfilebase);
first_ = all_(1);
cfg.yearbase = str2num(rawfilebase( first_-4:first_-1));

cfg = fillstruct(cfg, varargin, 0);


if (~strcmp(cfg.timebase, 'gga') & ~strcmp(cfg.timebase, 'log'))
   error('timebase must be ''log'' or ''gga''')
end

fext = 'gps'; %get info from gps files

outfile =     fullfile(cfg.outfpath, [rawfilebase, '.tim.gbin']);

if (exist(outfile, 'file') & (~cfg.redo))
   fprintf('not writing %s\n', outfile)
   return
end
  

fullrawfile = fullfile(cfg.rawfpath, [rawfilebase, '.raw']);
if (~exist(fullrawfile))
   error(sprintf('cannot find %s\n', fullrawfile))
end

%read the raw data
fprintf('reading raw data file %s/%s.raw ...\n', cfg.rawfpath, rawfilebase)
tic
data = read(cfg.instclass, fullrawfile, 'yearbase', cfg.yearbase,...
            'back_to_front', cfg.back_to_front);
toc

%read the log file and fix rawtimes to logger times
logfile = [fullrawfile, '.log']; 
logdata = textread(logfile); %measured unix times and offsets
meas_logger = logdata(:,1);
offsets = logdata(:,2);

if (meas_logger(1) > 10^6) %convert unix seconds to dday
   meas_logger = (meas_logger + to_sec(1970, 1, 1, 0, 0, 0) - ...
                               to_sec(cfg.yearbase, 1, 1, 0, 0, 0))/86400;
end
%this will be the same size as data.dday but with fixed meas_logger values
logger_dday = NaN*data.dday;   %need this many; but no more
if (length(data.dday) < length(meas_logger))     %bad checksum: lost data end
   meas_logger = meas_logger(1:length(data.dday));
   offsets = offsets(1:length(data.dday));
end
%index into nb dday, 
nominal_offset = median(diff(offsets));
iinto_dday = (offsets/nominal_offset) + 1;
%except there may be problems with the offsets
badi =  find(   (abs(iinto_dday - floor(iinto_dday)) ~= 0) | ...
                iinto_dday > length(data.dday));
iinto_dday(badi)=[];
meas_logger(badi)=[];

%now only use days from data.dday that have logger_dday 
logger_dday(iinto_dday) = meas_logger;
coef = polyfit(data.dday(iinto_dday),  logger_dday(iinto_dday), 1);
loggercorrected_rawdday = polyval(coef, data.dday);
% smoothr instrument times: with a 5-pt boxfilt on the time differences
dt = filt242(diff(loggercorrected_rawdday));
lc_rawddaymean = mean(loggercorrected_rawdday);  
filt_loggercorrected_rawdday = loggercorrected_rawdday(1) + cumsum([0; dt]);
%this might only be necessary for small files
flc_rawdday = filt_loggercorrected_rawdday - ...
    mean(filt_loggercorrected_rawdday) + lc_rawddaymean;
flc_rawdday = unscootch(flc_rawdday, loggercorrected_rawdday, cfg.scootch_sec);


%so now we have 
% flc_rawdday     = best version of raw dday corrected to logger times
% logger_dday     = logger times for flc_rawdday
% need to match that to gga (or the reverse)

%get the logger time info from the gps files
logger_range = [flc_rawdday(1) flc_rawdday(end)];
%bindata is [logger_dday, fullggadday, lon, lat]
bindata = bin_extract(cfg.binfpath, fext,  logger_range);

if (min(diff(flc_rawdday)) <= 0)
   warning('nonmonotonic filtered logger_corrected adcp times')
   fprintf(fid,...
           'nonmonotonic filtered logger_corrected adcp times in %s\n',...
           rawfilebase);
   keyboard
end

if (min(diff(flc_rawdday)) <= 0)
   warning('nonmonotonic logger times from binfile')
   keyboard
end

ggadday = table1([bindata(1,:)' bindata(2,:)'], flc_rawdday);
if (cfg.timebase == 'log')
   bestdday = flc_rawdday;   %use the filtered, logger_corrected adcp ddays
else %table from the filtered, logger_corrected adcp ddays to gga ddays
   bestdday = ggadday;
end
% logger time - gga
sec_pc_minus_utc = (flc_rawdday - ggadday)*86400;

gpsdata =  [flc_rawdday   bestdday sec_pc_minus_utc]';


if (length(gpsdata(1,:)) ~= length(data.dday))
   warning('mismatched gps times and data.dday')
   keyboard
end


%write the file
write_bin(outfile,gpsdata,'gps',...
          {'logger_dday', 'bestdday', 'sec_pc_minus_utc'})

if (cfg.verbose)
   fprintf('logger_range : %s - %s\n',...
           to_date(2001,logger_range(1)),...
           to_date(2001,logger_range(2)))
   subplot(211)
   plot(86400*diff(flc_rawdday))
   title('filtered logger_corrected adcp times')   
   yl(0,4)
   subplot(212)
   plot(86400*diff(table1([bindata(1,:)' bindata(2,:)'], flc_rawdday)))
   title('gga times if used')
   yl(0,4)
   drawnow
end

fid = fopen('gps2gbin.errs','a');
ymdhms = clock;
fprintf(fid,'writing %s at %s\n',  outfile, ...
        to_date(ymdhms(1), to_day(ymdhms(1), ymdhms),'s'));


if (cfg.timebase ~= 'log') & (min(diff(bestdday)) <= 0)
   warning(sprintf('nonmonotonic best (%s) adcp times', cfg.timebase))
   paused(2)
end


fclose(fid);

%%%%%%%%
function unscootched = unscootch(scootched, orig_dday, scootch_sec);
%
unscootched = scootched;

diffsec = 86400*(scootched - orig_dday);
scootchi = find( (abs(diffsec) > scootch_sec) & (diffsec > scootch_sec));
if (~isempty(scootchi))
   unscootched(scootchi) = orig_dday(scootchi) + scootch_sec/86400;
end
scootchi = find( (abs(diffsec) > scootch_sec) & (diffsec < - scootch_sec));
if (~isempty(scootchi))
   unscootched(scootchi) = orig_dday(scootchi) - scootch_sec/86400;
end


