%% main function
function nc_loaddap_odv_prof
urlbase = 'http://dods.mbari.org/cgi-bin/nph-nc/data/ssdsdata/deployments/m0/';
instrument_name = {'CTD0m','ts'};%,'ADCP','ecoflntu'x};%'imCTD';'ADCP-QC';'ecofl1';'ecofl2';'ecobb';'ocr3000';'gps';;;},...

addpath('/usr/local/bin');

[cruise, data,  header, time] = nc_loaddap (urlbase, instrument_name);
%"-15" prevents the last few hours of data from being written because they often are NaNs as a
%result of interpolation
datalines = length(data)-15;
data = data(1:datalines,:);
cruise = cruise(1:datalines,:);

fidcfg_old = fopen('/ssdsdata/deployments/m0/ascii_all_dep/M0PROF.cfg', 'r');
if fidcfg_old ~= -1,
    A = fread(fidcfg_old, '*char');
    fclose(fidcfg_old);
    A = A';
    A = A(3:end);
    prev_datalines = str2double(A);
    %only append the new data to the text files.  
    new_data = data(prev_datalines + 1 : datalines,:);
    new_cruise = cruise(prev_datalines + 1 : datalines,:);
    data = new_data;
    cruise = new_cruise;
end
    
if isempty(data),
    return,
end

%fidcfg = fopen('C:\Documents and Settings\sethb\My Documents\BOG\matlab\M0PROF.cfg', 'wt');
fidcfg = fopen('/ssdsdata/deployments/m0/ascii_all_dep/M0PROF.cfg', 'wt');
fprintf(fidcfg, '%s', ['//' int2str(datalines)]);
fclose(fidcfg);

if fidcfg_old ~= -1,
    fid = fopen('/ssdsdata/deployments/m0/ascii_all_dep/M0PROF.txt', 'r+');
    %fid = fopen('C:\Documents and Settings\sethb\My Documents\BOG\matlab\M0PROF.txt', 'r+');
    fprintf(fid, '%s', ['//' int2str(datalines)])
else
    fid = fopen('/ssdsdata/deployments/m0/ascii_all_dep/M0PROF.txt', 'wt');
    fprintf(fid, '%s\n', ['//' int2str(datalines)]);
    for h = 1:length(header),
        if h~= length(header),
            fprintf(fid, '%s\t', char(header{h}));
        else
            fprintf(fid, '%s\n', char(header{h}));
        end
    end
end
fclose(fid);

%dlmwrite('C:\Documents and Settings\sethb\My Documents\BOG\matlab\M0SURFtest2.txt', data, 'delimiter','\t','-append', 'roffset', 2,'precision','%.4f');
tblwrite_smb(data, header, cruise, '/ssdsdata/deployments/m0/ascii_all_dep/M0PROF.txt','\t',4,1)
%tblwrite_smb(data, header, cruise,'C:\Documents and Settings\sethb\My Documents\BOG\matlab\M0PROF.txt','\t',4,1)

end

%% Data Arrays

function [cruise, data,  header, time] = nc_loaddap (urlbase, instrument_name)
clear global time_hourly data_hourly
global time_hourly data_hourly

max_num_depths = 1;
depths = 0;

for i = 1:length(instrument_name),
    inst = char(instrument_name(i));

    %using the name of the instrument (inst) and the url base, load_data
    %retrieves the variable names and data contained in the netcdf files.
    %Data is obtained using loaddap and returned in a structure array
    %(inst_data).
    [variable_names, inst_data] = load_data (inst, urlbase);
    
    %loops through the multiple netcdf files for each instrument
     for f = 1:length(inst_data.(inst)),
        
        %loops through each field name for each netcdf file
        for v = 1:length(variable_names.(inst){f}),
            
            attribute = char(variable_names.(inst){f}(v));

            %checks that the attribute returns a structure, this indicates
            %that it is a multiple dimension variable and contains data
            attribute_struct = inst_data.(inst){f}.(attribute);
            if ~isstruct(attribute_struct),
                continue,
            end
            
            %puts the different fields within the attribute structure into
            %a cell array.  The first contains data, the second contains
            %time markers
            attribute_names = fieldnames(attribute_struct);
            
            %reads data into a cell/structure array for each attribute.
            %concatenates the data from each attribute across netcdf files.
            %Checks the size to determine if there are multiple data points
            %for each time, this indicates multiple depths.  
            if f==1,
                time_esecsx{i}.(attribute) = inst_data.(inst){f}.(attribute).(char(attribute_names(2)));
                temp = squeeze(inst_data.(inst){f}.(attribute).(char(attribute_names(1))));

                [m, n] = size(temp);

                if n == length(inst_data.(inst){f}.(attribute).(char(attribute_names(2)))),
                    temp = temp';
                end
                out_data{i}.(attribute) = temp;
                
                depth.(inst){f}.(attribute) = inst_data.(inst){f}.(attribute).(char(attribute_names(3)));
                if length(depth.(inst){f}.(attribute))>max_num_depths,
                    max_num_depths = length(depth.(inst){f}.(attribute));
                    depths = depth.(inst){f}.(attribute);
                end
            else
                time_esecsx{i}.(attribute) = [time_esecsx{i}.(attribute); inst_data.(inst){f}.(attribute).(char(attribute_names(2)))];
                temp = squeeze(inst_data.(inst){f}.(attribute).(char(attribute_names(1))));

                [m, n] = size(temp);

                if n == length(inst_data.(inst){f}.(attribute).(char(attribute_names(2)))),
                    temp = temp';
                end
                out_data{i}.(attribute) = [out_data{i}.(attribute); temp];
                
               depth.(inst){f}.(attribute) = inst_data.(inst){f}.(attribute).(char(attribute_names(3)));
                if length(depth.(inst){f}.(attribute))>max_num_depths,
                    max_num_depths = length(depth.(inst){f}.(attribute));
                    depths = depth.(inst){f}.(attribute);
                end

            end

        end
    end
    
    %interpolates data to an hourly grid
    time_grid (time_esecsx, out_data, i, inst);
        
end
% 
% for v = 1:length(instrument_name),
%     variable_names = fieldnames(data_hourly.(instrument_name{v}));
%     for n = 1:length(variable_names),
%         temp = data_hourly.(instrument_name{v}).(variable_names{n});
%         [d, e] = size(temp);
%         if e<length(depths),
%            temp(:,e+1:length(depths)) = NaN;
%         end
%         data_hourly.(instrument_name{v}).(variable_names{n}) = temp;
%             
% %             for j = 1:length(depths),
% %                 match = find(depth.(instrument_name{v}){1}.(variable_names{n})(j)==depths);
% %                 
% %                 data_hourly.(instrument_name{v}).(variable_names{n})(:,j) = data_hourly.(instrument_name{v}).(variable_names{n})(:,match);
% %                
% %             end
%         end
%         
%         
%     end
% end

time_edays = time_hourly/(3600*24);
epoch_start = datenum('1-jan-1970');
time = epoch_start + time_edays;

mdy = datestr(time, 'mm/dd/yyyy');
hhmm = datestr(time, 'HH:MM');

lon = inst_data.(inst){f}.(attribute).(char(attribute_names(5)));
lat = inst_data.(inst){f}.(attribute).(char(attribute_names(4)));
%sets "Bot.Depth" to 0, purely for file structure benefit
%depth(1:length(time),1) = 0;
days_2003 = time - (datenum('1-Jan-2003')) + 1;


for l = 1:length(instrument_name),
    data_names = fieldnames(data_hourly.(char(instrument_name{l})));
    for m = 1:length(data_names),
        miss_val = find(data_hourly.(char(instrument_name{l})).(char(data_names{m}))(:,:)<-100);
        data_hourly.(char(instrument_name{l})).(char(data_names{m}))(miss_val)=NaN;
    end
end


for row = 1:length(time),
    for dep = 1:length(depths),
        cruise{(row-1)*length(depths)+dep,1} = 'M0';
        cruise{(row-1)*length(depths)+dep,2} = num2str(row); %num2str((row-1)*length(depths)+dep);
        cruise{(row-1)*length(depths)+dep,3} = '*';
        cruise{(row-1)*length(depths)+dep,4} = mdy(row,1:10);
        cruise{(row-1)*length(depths)+dep,5} = hhmm(row,1:5);

        data((row-1)*length(depths)+dep,1) = lon+360;
        data((row-1)*length(depths)+dep,2) = lat;
        data((row-1)*length(depths)+dep,3) = depths(dep);
        data((row-1)*length(depths)+dep,4) = days_2003(row);
        data((row-1)*length(depths)+dep,5) = depths(dep); %for now sensor depth is all 0, since this is the surface file...
        data((row-1)*length(depths)+dep,6) = 0;
        
        if dep == 1,
            data((row-1)*length(depths)+dep,7) = data_hourly.CTD0m.Temperature(row,dep);
            data((row-1)*length(depths)+dep,9) = data_hourly.CTD0m.Salinity(row,dep);
        else
            data((row-1)*length(depths)+dep,7) = data_hourly.ts.SEA_WATER_TEMPERATURE(row,dep);
            data((row-1)*length(depths)+dep,9) = data_hourly.ts.SEA_WATER_SALINITY(row,dep);
        end
        data((row-1)*length(depths)+dep,8) = 0;
        data((row-1)*length(depths)+dep,10) = 0;
    end
end




header = {'Cruise', 'Station', 'Type', 'mon/day/yr', 'hh:mm',...
    'Longitude [degrees_east]', 'Latitude [degrees_north]', 'Bot.Depth [m]',...
    'DaysSince2003[d]','SensorDepth[m]','QF', 'Temperature[�C]', 'QF', 'Salinity', 'QF'};
    
end

%% Load Variables
function [variable_names, instrument] = load_data(inst, urlbase)

switch char(inst);

    case {'CTD0m'}
        file.(inst){1} = '200405/m0_ctd0000_20040727_original.nc';
        file.(inst){2} = '200410/m0_ctd0000_20041006_original.nc';
        file.(inst){3} = '200504/m0_ctd0001_20050418_original.nc';
        file.(inst){4} = '200607/m0_ctd0001_20060731_original.nc';
        file.(inst){5} = '200701/m0_ctd0000_20070130_original.nc';
        file.(inst){6} = '200706/m0_ctd0001_20070621_original.nc';
    case {'imCTD'}
        file.(inst){1} = '200405/m0_imctd_20040526_original.nc';
        file.(inst){2} = '200410/m0_imctd_20041006_original.nc';
        file.(inst){3} = '200504/m0_imctd_20050418_original.nc';
        file.(inst){4} = '200607/m0_imctd_20060727_original.nc';
        file.(inst){5} = '200701/m0_imctd_20070130_original.nc';
        file.(inst){6} = '200706/m0_imctd_20070621_original.nc';
    case{'hourly'}
        file.(inst){1} = '200405/hourlyM0_20040525.nc';
        file.(inst){2} = '200410/hourlyM0_20041006.nc';
        file.(inst){3} = '200504/hourlyM0_20050418.nc';
        file.(inst){4} = '200607/hourlyM0_20060731.nc';
        file.(inst){5} = '200701/hourlyM0_20070130.nc';
        file.(inst){6} = '200706/hourlyM0_.nc';
    case{'metsys'}
        file.(inst){1} = '200405/m0_metsys_20040525_original.nc';
        file.(inst){2} = '200410/m0_metsys_20041006_original.nc';
        file.(inst){3} = '200504/m0_metsys_20050418_original.nc';
        file.(inst){4} = '200607/m0_metsys_20060801_original.nc';
        file.(inst){5} = '200607/m0_metsys_20061108_original.nc';
        file.(inst){6} = '200701/m0_metsys_20070130_original.nc';
        file.(inst){7} = '200706/m0_metsys_20070621_original.nc';
    case{'CO2QC'}
        file.(inst){1} = 'qc/m0_pco2_qc.nc';
        file.(inst){2} = '200701/m0_pco2_20070130.nc';
        %mike mccann hasn't updated yet...
    case{'CO2'}
        file.(inst){1} = '200607/m0_pco2_20060731.nc';
        file.(inst){2} = '200701/m0_pco2_20070130.nc';
    case{'ADCP'}
        file.(inst){1} = '200405/m0_adcp1267_20040525.nc';
        file.(inst){2} = '200410/m0_adcp1267_20041006.nc';
        file.(inst){3} = '200504/m0_adcp1406_20050418.nc';
        file.(inst){4} = '200607/m0_adcp1267_20060731.nc';
        file.(inst){5} = '200701/m0_adcp1267_20070130.nc';
    case{'ecoflntu'}
        file.(inst){1} = '200607/m0_ecoflntu_20060731.nc';
        file.(inst){2} = '200607/m0_ecoflntu_20061108.nc';
        file.(inst){3} = '200701/m0_ecoflntu_20070130.nc';
        file.(inst){4} = '200706/m0_ecoflntu_20070621.nc';
    case{'ts'}
        file.(inst){1} = '200607/m0_ts_20060731.nc';
        file.(inst){2} = '200701/m0_ts_20070130.nc';
        file.(inst){3} = '200706/m0_ts_20070621.nc';
    case{'ocr3000'}
        file.(inst){1} = '200405/m0_ocr3000_prev.nc';
        file.(inst){2} = '200607/m0_ocr3000_20060721.nc';
        file.(inst){3} = '200701/m0_ocr3000_20070130.nc';
        file.(inst){4} = '200706/m0_ocr3000_20070621.nc';
    case{'ecofl1'}
        file.(inst){1} = '200405/m0_ecofl_20040526.nc';
    case{'ecofl2'}
        file.(inst){1} = '200410/m0_ecofl_20041006.nc';
        file.(inst){2} = '200504/m0_ecofl_20050418.nc';
    case{'ecobb'}
        file.(inst){1} = '200504/m0_ecobb_20050418.nc';
        file.(inst){2} = '200504/m0_ecobb_20050728.nc';
    case{'gps'}
        file.(inst){1} = '200405/m0_gps1313_20040525.nc';
        file.(inst){2} = '200410/m0_gps1313_20041006.nc';
        file.(inst){3} = '200504/m0_gps1313_20050418.nc';
        file.(inst){4} = '200607/m0_gps1564_20060731.nc';
        file.(inst){5} = '200701/m0_gps1564_20070130.nc';
        file.(inst){6} = '200706/m0_gps1313_20070621.nc';
    case{'ADCP-QC'}
        file.(inst){1} = 'qc/200406/m0_adcp_20040603_workhorse.nc';
    case{'metsys-QC'}
        file.(inst){1} = 'qc/200406/m0_metsys_20040603.nc';
    case{'CTD-QC'}
        file.(inst){1} = 'qc/200406/m0_tstring_20040603.nc';
    case{'chl'}
        file.(inst){1} = '../../ocr3000cimt/netcdf/200701/chl_490_Kd_0_10_mean.nc';
        file.(inst){2} = '../../ocr3000cimt/netcdf/200706/chl_490_Kd_0_10_mean.nc';
end

for f = 1:length(file.(inst)),
    variable_names.(inst){f} = fieldnames(loaddap([urlbase char(file.(inst){f})]));
    instrument.(inst){f} = loaddap([urlbase char(file.(inst){f})]);
end
end

%% time gridding

%retrieves the names of different data types that are now stored in the
%out_data cell/structure array

function  time_grid (time_esecsx, out_data, i, inst)
data_names = fieldnames(out_data{i});

global time_hourly data_hourly

for d = 1:length(data_names),
    time_esecs_start = time_esecsx{i}.(char(data_names(d)))(1);
    time_esecs_end = time_esecsx{i}.(char(data_names(d)))(end);

    %moves the start and end time to an even hour
    time_esecs_start = double((int32(time_esecs_start/3600)).*3600);
    time_esecs_end = double((int32(time_esecs_end/3600)+1).*3600);%the plus one establishes the end time after the last timestamp
    
    %sets the interval to which the data will be interpolated, in this
    %case 1 hour (60 min x 60 secs/min)
    interval = 3600;

    %sets the hourly grid to the start and end points of the first variable
    %it encounters.  Perhaps not the best plan, but a starting point.
    temp = time_esecs_start:interval:time_esecs_end;
    if isempty(time_hourly),
        time_hourly = temp;
    end
    
    %for some reason the time stamps aren't monotonic- possibly due to
    %overlap in deployments? this fixes that.
    y = diff([0;time_esecsx{i}.(char(data_names(d)))]);
    w = find(y>=0);
    mono_time = time_esecsx{i}.(char(data_names(d)))(w);
    mono_data = out_data{i}.(char(data_names(d)))(w,:);

    neg = 1;
    while (isempty(neg)==0),
        q = diff([0;mono_time]);
        w = find(q>0);
        neg = find(q<=0);
        mono_time = mono_time(w);
        mono_data = mono_data(w,:);
    end

    data_hourly.(inst).(char(data_names(d))) = interp1(mono_time, mono_data, time_hourly);
    
    [m, n, o] = size (data_hourly.(inst).(char(data_names(d))));

    if n==length(time_hourly),
        data_hourly.(inst).(char(data_names(d))) = (data_hourly.(inst).(char(data_names(d))))';
    end
    
end
end
