write_FVCOM_elevtide

 Write an FVCOM surface elevation time series forcing file 

 write_FVCOM_elevtide(Mobj,MJD,ElevationFile,MyTitle)

 DESCRIPTION:
    Write an FVCOM NetCDF surface elevation forcing file

 INPUT:
   Mobj         = Matlab mesh object.
   MJD          = list of modified Modified Julian Dates of size [times]
                   (defined as unlimited in the NetCDF file).
   ElevationFile    = name of NetCDF file.
   MyTitle      = casename title, written as global attribute of NetCDF file.

 OUTPUT:
    ElevationFile, A NetCDF FVCOM surface elevations tide forcing file

 EXAMPLE USAGE
    write_FVCOM_elevtide(Mobj,MJD,ElevationFile,MyTitle)

 Author(s):  
    Pierre Cazenave (Plymouth Marine Laboratory)
    Karen Thurston (National Oceanography Centre Liverpool)
 
 Revision history
    2012-08-08 (PWC) First version.
    2012-11-14 (PWC) Updated to expect Modified Julian Day rather than doing 
    the conversion in here. Also put the pieces in set_elevtide in here to
    simplify the process of writing out an elevation input file.
    2012-12-04 (KJT) Updated to use surface elevation and open boundary 
    nodes from Mobj.
   
==============================================================================

0001 function write_FVCOM_elevtide(Mobj,MJD,ElevationFile,MyTitle)
0002 % Write an FVCOM surface elevation time series forcing file
0003 %
0004 % write_FVCOM_elevtide(Mobj,MJD,ElevationFile,MyTitle)
0005 %
0006 % DESCRIPTION:
0007 %    Write an FVCOM NetCDF surface elevation forcing file
0008 %
0009 % INPUT:
0010 %   Mobj         = Matlab mesh object.
0011 %   MJD          = list of modified Modified Julian Dates of size [times]
0012 %                   (defined as unlimited in the NetCDF file).
0013 %   ElevationFile    = name of NetCDF file.
0014 %   MyTitle      = casename title, written as global attribute of NetCDF file.
0015 %
0016 % OUTPUT:
0017 %    ElevationFile, A NetCDF FVCOM surface elevations tide forcing file
0018 %
0019 % EXAMPLE USAGE
0020 %    write_FVCOM_elevtide(Mobj,MJD,ElevationFile,MyTitle)
0021 %
0022 % Author(s):
0023 %    Pierre Cazenave (Plymouth Marine Laboratory)
0024 %    Karen Thurston (National Oceanography Centre Liverpool)
0025 %
0026 % Revision history
0027 %    2012-08-08 (PWC) First version.
0028 %    2012-11-14 (PWC) Updated to expect Modified Julian Day rather than doing
0029 %    the conversion in here. Also put the pieces in set_elevtide in here to
0030 %    simplify the process of writing out an elevation input file.
0031 %    2012-12-04 (KJT) Updated to use surface elevation and open boundary
0032 %    nodes from Mobj.
0033 %
0034 %==============================================================================
0035 
0036 global ftbverbose 
0037 report = false;
0038 if(ftbverbose); report = true; end
0039 subname = 'write_FVCOM_elevtide';
0040 if(report); fprintf('\n'); end
0041 if(report); fprintf(['begin : ' subname '\n']); end
0042 
0043 % Get a list of the open boundary nodes. Transpose Mobj.obc_nodes so the
0044 % order of the boundary nodes is preserved.
0045 tmpObcNodes = Mobj.obc_nodes';
0046 % Flip it back so it's the same shape as it would have been using the old
0047 % code.
0048 ObcNodes = tmpObcNodes(tmpObcNodes~=0)';
0049 
0050 %------------------------------------------------------------------------------
0051 % Sanity check on input and dimensions
0052 %------------------------------------------------------------------------------
0053 nTimes = numel(MJD);
0054 if(report); fprintf('Number of time steps %d\n',nTimes); end
0055 
0056 nObcs = numel(ObcNodes);
0057 if(report); fprintf('Number of Open Boundary Nodes %d\n',nObcs); end
0058 
0059 [chk1, chk2] = size(Mobj.surfaceElevation);
0060 if nObcs ~= chk1 || nTimes ~= chk2
0061     fprintf('Surface elevation dimensions do not match time series and number of boundary nodes.\n')
0062     fprintf('Surface elevation nodes and time sizes: (%d, %d)\n', chk1, chk2)
0063     fprintf('Boundary nodes size: %d\n', nObcs)
0064     fprintf('Times size: %d\n', nTimes)
0065     error('Input data sizes do not match. Check and try again.');
0066 end
0067 
0068 %%
0069 %------------------------------------------------------------------------------
0070 % Dump the file
0071 %------------------------------------------------------------------------------
0072 
0073 nc=netcdf.create(ElevationFile,'clobber');
0074 
0075 % define global attributes
0076 netcdf.putAtt(nc,netcdf.getConstant('NC_GLOBAL'),'type','FVCOM TIME SERIES ELEVATION FORCING FILE')
0077 netcdf.putAtt(nc,netcdf.getConstant('NC_GLOBAL'),'title',MyTitle)
0078 netcdf.putAtt(nc,netcdf.getConstant('NC_GLOBAL'),'history','FILE CREATED using write_FVCOM_elevtide')
0079 
0080 % define dimensions
0081 nobc_dimid=netcdf.defDim(nc,'nobc',nObcs);
0082 time_dimid=netcdf.defDim(nc,'time',netcdf.getConstant('NC_UNLIMITED'));
0083 date_str_len_dimid=netcdf.defDim(nc,'DateStrLen',26);
0084 
0085 % define variables and attributes
0086 nobc_varid=netcdf.defVar(nc,'obc_nodes','NC_INT',nobc_dimid);
0087 netcdf.putAtt(nc,nobc_varid,'long_name','Open Boundary Node Number');
0088 netcdf.putAtt(nc,nobc_varid,'grid','obc_grid');
0089 
0090 iint_varid=netcdf.defVar(nc,'iint','NC_INT',time_dimid);
0091 netcdf.putAtt(nc,iint_varid,'long_name','internal mode iteration number');
0092 
0093 time_varid=netcdf.defVar(nc,'time','NC_FLOAT',time_dimid);
0094 netcdf.putAtt(nc,time_varid,'long_name','time');
0095 netcdf.putAtt(nc,time_varid,'units','days since 1858-11-17 00:00:00');
0096 netcdf.putAtt(nc,time_varid,'format','modified julian day (MJD)');
0097 netcdf.putAtt(nc,time_varid,'time_zone','UTC');
0098 
0099 itime_varid=netcdf.defVar(nc,'Itime','NC_INT',time_dimid);
0100 netcdf.putAtt(nc,itime_varid,'units','days since 1858-11-17 00:00:00');
0101 netcdf.putAtt(nc,itime_varid,'format','modified julian day (MJD)');
0102 netcdf.putAtt(nc,itime_varid,'time_zone','UTC');
0103 
0104 itime2_varid=netcdf.defVar(nc,'Itime2','NC_INT',time_dimid);
0105 netcdf.putAtt(nc,itime2_varid,'units','msec since 00:00:00');
0106 netcdf.putAtt(nc,itime2_varid,'time_zone','UTC');
0107 
0108 Times_varid=netcdf.defVar(nc,'Times','NC_CHAR',[date_str_len_dimid, time_dimid]);
0109 netcdf.putAtt(nc,Times_varid,'time_zone','UTC');
0110 
0111 elevation_varid=netcdf.defVar(nc,'elevation','NC_FLOAT',[nobc_dimid, time_dimid]);
0112 netcdf.putAtt(nc,elevation_varid,'long_name','Open Boundary Elevation');
0113 netcdf.putAtt(nc,elevation_varid,'units','meters');
0114 
0115 % end definitions
0116 netcdf.endDef(nc);
0117 
0118 % write data
0119 netcdf.putVar(nc,nobc_varid,ObcNodes);
0120 netcdf.putVar(nc,iint_varid,0,nTimes,1:nTimes);
0121 netcdf.putVar(nc,time_varid,0,nTimes,MJD);
0122 netcdf.putVar(nc,itime_varid,floor(MJD));
0123 netcdf.putVar(nc,itime2_varid,0,nTimes,mod(MJD,1)*24*3600*1000);
0124 nStringOut = char();
0125 for i=1:nTimes
0126     [nYr, nMon, nDay, nHour, nMin, nSec] = mjulian2greg(MJD(i));
0127     if strcmp(sprintf('%02i', nSec), '60')
0128         % Fix some weirdness with mjulian2greg. I think this is caused by
0129         % rounding errors. My testing suggests this is not a problem around
0130         % midnight, so the number of days (and thus possibly months and
0131         % years) is unaffected.
0132         if mod(nMin + 1, 60) == 0
0133             % Up the hour by one too
0134             nHour = mod(nHour + 1, 24);
0135         end
0136         nMin = mod(nMin + 1, 60);
0137         nSec = 0;
0138     end
0139     nDate = [nYr, nMon, nDay, nHour, nMin, nSec];
0140     nStringOut = [nStringOut, sprintf('%04i/%02i/%02i %02i:%02i:%02i       ',nDate)];
0141 end
0142 netcdf.putVar(nc,Times_varid,nStringOut);
0143 netcdf.putVar(nc,elevation_varid,Mobj.surfaceElevation);
0144 
0145 % close file
0146 netcdf.close(nc);
0147 
0148 if(report); fprintf(['end   : ' subname '\n']); end;
0149