function write_roms_otps_ncfile ( Tide, grid_ncfile, output_tides_ncfile,base_date,model_file)
% WRITE_ROMS_OTPS_NCFILE:
%
% PARAMETERS:
% Input:
%     Tide:
%        Structure with arrays to be written to file
%     grid_ncfile:
%        Existing grid file.  Will not be written to.
%     output_tides_ncfile:
%        To be created.
%    
%Revised by Eli Hunter 3/7/07 
%1) Include base_date from otps2frc_v2 input
%2) Include model_file from otps2frc_v2 input



fprintf ( 1, '%s:  writing %s...\n', mfilename, output_tides_ncfile );


%
% Create the output tidal file
grid = nc_info ( grid_ncfile );
eta_rho_info = nc_getdiminfo ( grid_ncfile, 'eta_rho' );
xi_rho_info = nc_getdiminfo ( grid_ncfile, 'xi_rho' );
lon_rho = nc_varget ( grid_ncfile, 'lon_rho' );
lat_rho = nc_varget ( grid_ncfile, 'lat_rho' );
mask_rho = nc_varget ( grid_ncfile, 'mask_rho' );


%
[ncid, status] = mexnc ( 'create', output_tides_ncfile, nc_clobber_mode );
if ( status ~= 0 )
	msg = sprintf ( '%s:  mexnc:create failed on %s.\n', mfilename, output_tides_ncfile );
	error ( msg );
end

[two_dimid, status] = mexnc ( 'DEF_DIM', ncid, 'two', 2 );
if ( status ~= 0 )
	msg = sprintf ( '%s:  %s.\n', mfilename, mexnc ( 'strerror', status ) );
	error ( msg );
end

status = mexnc ( 'close', ncid );
if ( status ~= 0 )
	msg = sprintf ( '%s:  mexnc:close failed on %s.\n', mfilename, output_tides_ncfile );
	error ( msg );
end

 nc_add_dimension ( output_tides_ncfile, 'eta_rho', eta_rho_info.Length );
 nc_add_dimension ( output_tides_ncfile, 'xi_rho', xi_rho_info.Length );
nc_add_dimension ( output_tides_ncfile, 'tide_period', 0 );


%
% Now define the variables
v.Name = 'lat_rho';
v.Dimension = { 'eta_rho', 'xi_rho' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'latitude of RHO-points';
v.Attribute(2).Name = 'units';
v.Attribute(2).Value = 'degree_north';
v.Attribute(3).Name = 'field';
v.Attribute(3).Value = 'lat_rho, scalar';
nc_addvar ( output_tides_ncfile, v );


v.Name = 'lon_rho';
v.Dimension = { 'eta_rho', 'xi_rho' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'longitude of RHO-points';
v.Attribute(2).Name = 'units';
v.Attribute(2).Value = 'degree_east';
v.Attribute(3).Name = 'field';
v.Attribute(3).Value = 'lon_rho, scalar';
nc_addvar ( output_tides_ncfile, v );


v.Name = 'mask_rho';
v.Dimension = { 'eta_rho', 'xi_rho' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'mask on RHO-points';
v.Attribute(2).Name = 'option_0';
v.Attribute(2).Value = 'land';
v.Attribute(3).Name = 'option_1';
v.Attribute(3).Value = 'water';
nc_addvar ( output_tides_ncfile, v );


v.Name = 'tide_period';
v.Dimension = { 'tide_period' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'tide angular period';
v.Attribute(2).Name = 'units';
v.Attribute(2).Value = 'hours';
v.Attribute(3).Name = 'field';
v.Attribute(3).Value = 'tide_period, scalar, series';

nc_addvar ( output_tides_ncfile, v );




v.Name = 'tide_Ephase';
v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'tidal elevation phase angle';
v.Attribute(2).Name = 'units';
v.Attribute(2).Value = 'degrees, time of maximum elevation with respect to chosen time origin';
v.Attribute(3).Name = 'field';
v.Attribute(3).Value = 'tide_Ephase, scalar, series';

nc_addvar ( output_tides_ncfile, v );





v.Name = 'tide_Eamp';
v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'tidal elevation amplitude';
v.Attribute(2).Name = 'units';
v.Attribute(2).Value = 'meter';
v.Attribute(3).Name = 'field';
v.Attribute(3).Value = 'tide_Eamp, scalar, series';

nc_addvar ( output_tides_ncfile, v );



	
	
	v.Name = 'tide_Cphase';
	v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
	v.Attribute(1).Name = 'long_name';
	v.Attribute(1).Value = 'tidal current phase angle';
	v.Attribute(2).Name = 'units';
	v.Attribute(2).Value = 'degrees, time of maximum velocity with respect chosen time origin';
	v.Attribute(3).Name = 'field';
	v.Attribute(3).Value = 'tide_Cphase, scalar';
	
	nc_addvar ( output_tides_ncfile, v );
	
	
	
	
	v.Name = 'tide_Cangle';
	v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
	v.Attribute(1).Name = 'long_name';
	v.Attribute(1).Value = 'tidal current inclination angle';
	v.Attribute(2).Name = 'units';
	v.Attribute(2).Value = 'degrees between semi-major axis and East' ;
	v.Attribute(3).Name = 'field';
	v.Attribute(3).Value = 'tide_Cangle, scalar';
	
	nc_addvar ( output_tides_ncfile, v );
	
	
	v.Name = 'tide_Cmin';
	v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
	v.Attribute(1).Name = 'long_name';
	v.Attribute(1).Value = 'minimum tidal current, ellipse semi-minor axis';
	v.Attribute(2).Name = 'units';
	v.Attribute(2).Value = 'meter second-1' ;
	v.Attribute(3).Name = 'field';
	v.Attribute(3).Value = 'tide_Cmin, scalar';
	
	nc_addvar ( output_tides_ncfile, v );
	
	
	
	
	v.Name = 'tide_Cmax';
	v.Dimension = { 'tide_period', 'eta_rho', 'xi_rho' };
	v.Attribute(1).Name = 'long_name';
	v.Attribute(1).Value = 'maximum tidal current, ellipse semi-major axis';
	v.Attribute(2).Name = 'units';
	v.Attribute(2).Value = 'meter second-1' ;
	v.Attribute(3).Name = 'field';
	v.Attribute(3).Value = 'tide_Cmax, scalar';
	
	nc_addvar ( output_tides_ncfile, v );
	


clear v;
v.Name = 'tidal_constituents';
v.Nctype = 'char';
v.Dimension = { 'tide_period', 'two' };
v.Attribute(1).Name = 'long_name';
v.Attribute(1).Value = 'Tidal Constituent Names';

nc_addvar ( output_tides_ncfile, v );
	

%
% Now write the unlimited dimension data
output_buffer.tide_period = Tide.period;
output_buffer.tide_Ephase = Tide.Ephase;
output_buffer.tide_Eamp = Tide.Eamp;
output_buffer.tide_Cphase = Tide.Cphase;
output_buffer.tide_Cangle = Tide.Cangle;
output_buffer.tide_Cmin = Tide.Cmin;
output_buffer.tide_Cmax = Tide.Cmax;

%
% Construct the names
for j = 1:length(Tide.names)
	output_buffer.tidal_constituents(j,:) = Tide.names{j}
end


nc_addnewrecs ( output_tides_ncfile, output_buffer, 'tide_period' );

dname=input('Enter domain namefor NetCDf attribute:');

%
% And finally add some comments
Attribute(1).Name = 'type';
Attribute(1).Value = 'ROMS Forcing File';
Attribute(2).Name = 'title';
Attribute(2).Value =['Forcing for ' dname '  domain'];
Attribute(3).Name = 'base_date';
Attribute(3).Value = ['days since ' datestr(base_date,31)];
Attribute(4).Name = 'grid_file';
Attribute(4).Value = grid_ncfile;
Attribute(5).Name = 'source';
Attribute(5).Value = 'OTPS';
Attribute(6).Name = 'source_url';
Attribute(6).Value = 'http://www.coas.oregonstate.edu/research/po/research/tide/region.html';
Attribute(7).Name = 'history';
Attribute(7).Value = sprintf ( '%s:  Created by %s with %s.\n', datestr(now), getenv('USER'), mfilename );
Attribute(8).Name = 'comment';
Attribute(8).Value = ['Inputs for OTPS executable "extract_HC" ' ...
                    'created with m-file roms2ll.m using the grid ' ...
                    'file as input.' model_file ' was used as the OTPS regional model.   '];

%
% Name the tidal constituents.
const_string = Tide.names{1};
for j = 2:length(Tide.names)
	const_string = sprintf ( '%s, %s', const_string, Tide.names{j} );
end
Attribute(9).Name = 'tidal_constituents';
Attribute(9).Value = const_string;

for j = 1:length(Attribute)
	nc_attput ( output_tides_ncfile, nc_global, Attribute(j).Name, Attribute(j).Value );
end



%
% transfer the lon and lats
lon = nc_varget ( grid_ncfile, 'lon_rho' );
nc_varput ( output_tides_ncfile, 'lon_rho', lon );
lat = nc_varget ( grid_ncfile, 'lat_rho' );
nc_varput ( output_tides_ncfile, 'lat_rho', lat );
mask_rho = nc_varget ( grid_ncfile, 'mask_rho' );
nc_varput ( output_tides_ncfile, 'mask_rho', mask_rho );

