% Script CDF_DEMO
% demonstrates creation and manipulation of netCDF time series 
% files with EPIC conventions
%   Uses create_ts.m to create netCDF file
%   Uses getts.m to extract information from the created file
%

%
% Rich Signell               |  rsignell@crusty.er.usgs.gov
% U.S. Geological Survey     |  (508) 457-2229  |  FAX (508) 457-2310
% Quissett Campus            | 
% Woods Hole, MA  02543      |
%
%%%%%%%%%%%%%%%%%%%%%
% PART A. CREATING THE NETCDF TIME SERIES DATA
%%%%%%%%%%%%%%%%%%%%%
%1.  File Name
  cdf_name= 'flux.cdf';   
 
%2.  Data matrix (each column a time series)
  load flux.dat
  var=flux;

%3.  Variable names and units
  var_names=['shortwave_radiation '
             'sensible_heat_flux  '
             'latent_heat_flux    '
             'longwave_radiation  '
             'total_heat_flux     '];

  units=['W/m2'
         'W/m2'
         'W/m2'
         'W/m2'
         'W/m2' ];

%3.  Lat,Lon,Depth of time series data
  lat=  41.525;
  lon=  -70.67333;
  idepth=0.;     % instrument depth in meters

%4.  Start time and time interval
  start= [1990 1 2 17 0 0];     %[yyyy mm dd hr mi sc]
  dt= 3600*24;    % delta t in seconds

%5.  Water Depth 
  wdepth=30.;   % water depth in meters

%6.  Comment
  cmt=[
  'Latent and Sensible heat flux calculated by Bill Large routine bulk.f  ' 
  'using relative humidity from Logan and wind speed, air temp, sst and   '
  'air pressure from Boston Buoy B (NODC 44013).                          '];
%
%7. Call make_ts

make_ts(cdf_name,var,var_names,lat,lon,start,dt,wdepth,idepth,units,cmt');

%%%%%%%%%%%%%%%%%%%%%
% PART B. READING THE NETCDF TIME SERIES DATA
%%%%%%%%%%%%%%%%%%%%%

% 1. Use mcvgt to get all the sensible heat flux. 
%    mcvgt will work with 1D data, but also with 2D, 3D or arbitrary 
%    dimensioned data.  It doesn't require any particular conventions.  
%    You can specify the "corner" and "edges" of the "hyperslab" (netCDF
%    lingo) if you don't want all the data.

var='sensible_heat_flux';
disp(['read ' var ' with mcvgt and plot it']);
sens=mcvgt(cdf_name,var);   % get data from netCDF file
plot(sens)                                   % plot it
title(var)
disp(' hit [return] to continue')
pause

% 2. Use getts to get all shortwave radiation, time vector and long_name

var='shortwave_radiation';
disp(['read ' var ' with getts and plot vs. Gregorian time']);
[sw,jd,long_name]=getts(cdf_name,var); 
jd0=julian([1990 1 1 0 0 0]);
plot(jd-jd0,sw)       
gregaxis
title(long_name)
disp(' hit [return] to continue')
pause

% 3. Use getts to get latent heat flux between Jan 2, 1991 and Feb 28, 1991.
%    Get auxiliary information (lat,lon,depth, units, etc.) and use for plot
var='latent_heat_flux';
disp(['read ' var ' from Jan 2 to Feb 28, 1991 and plot vs. Gregorian Time']);
[latent,jd,long_name,units,lat,lon,idepth,wdepth]=...
  getts(cdf_name,'latent_heat_flux',[1991 1 2 0 0 0],[1991 2 28 0 0 0]);
plot(jd-jd0,latent)     
gregaxis
title([long_name ' at ' num2str(lon) ',' num2str(lat)])
ylabel(units)
