HOPS PE model initialization and update PE_INITIAL package: HOST: maelstrom.harvard.edu (128.103.2.50) TAR FILES: pub/HOPS/PE_initial/peini_9.10.tar.Z (252955:1034752 bytes) pub/HOPS/PE_initial/Ex_peini_9.10.tar.gz (10259707:20267520 bytes) VERSION: 9.10 (August 17, 2001) ORIGIN: Harvard University, Cambridge Massachusetts Harvard Ocean Prediction System (HOPS) DEVELOPERS: Patrick J. Haley Jr. (haley@pacific.harvard.edu) Carlos J. Lozano (lozano@pacific.harvard.edu) LIBRARIES: (1) NetCDF, version 3.3.1 University Corporation for Atmospheric Research/UNIDATA http://www.unidata.ucar.edu/ (2) GNUmake version 3.75 Free Software Foundation http://www.gnu.ai.mit.edu/ ============ INTRODUCTION ============ This directory contains HOPS's Primitive Equation (PE) model initialization and update package. This package is used to prepare fields for use in the PE model. Its primary tasks are: - Interpolating fields in the vertical to the PE grid. - Constructing estimates of the internal mode velocity (either from dynamic height, streamfunction or thermal wind relation) - Constructing estimates of the transport streamfunction (according to user-specified instructions) The PE model requires at initialization 3-D fields of density and velocity. The density is provided in terms of temperature (degree Celsius) and salinity (PSU) whereas the velocity is given in terms of internal (baroclinic) and external (barotropic) components (cm/s). The external mode component is given in terms of a transport (vertically integrated) streamfunction (cm^3/s). If applicable, it provides the additional tracers needed to run the biological module. This package also provides the boundary conditions for the above fields. In addition, this package provides the bottom topography, Land/Sea mask data (if any), mean T-S profile used in the computation of the mean density, and the necessary parameters to set-up the model horizontal and vertical grids. The initial, boundary, and update (assimilation) fields are written into a NetCDF file which is read by the PE model. A second executable, make_transbdy, is included to create boundary conditions for constructing the transport streamfunction by interpolation from the results for a larger domain. This is primarily motivated by the desire to keep the initialization fields for nested domains in agreement. ============ Installation ============ This package is available over the INTERNET via anonymous FTP from maelstrom.harvard.edu (128.103.2.50). When connected, you will be in the FTP directory. To obtain this package, go to the directory "pub/HOPS/PE_initial" and get files: Readme.peini This file. peini_9.10.tar.Z Compressed tar file of the PE_INITIAL package. Ex_peini_9.10.tar.gz Compressed tar file of an example using the PE_INITIAL package. To install this package, simply go to the directory in which you want to put the PE_INITIAL package and execute the following commands: zcat peini_9.10.tar.Z | tar -pxvf - gzip -dc Ex_peini_9.10.tar.gz | tar -pxvf - The tar file peini_9.10.tar.Z contains the following files: PE_initial.F make_transbdy.F all_uc.F avg_isle.F bess1d.F bess2d.F bess2d_pmsk.F blank.F blkdat.F cgpois.F cgpois1.F check_depths.F day_code.F defcdf.F exitus.F external.F fdiagn.F filter.F fsave.F get_cdfdat.F get_data.F get_date.F get_dom.F get_ewpt.F get_isval.F get_nbhrs.F get_opn.F get_parm.F get_pbar.F get_trc.F get_unif.F head_unif.F header_fm.F header_qg.F headln.F inside.F interp.F intrp2d.F intrp_pbar.F isincrs.F length.F lintrp.F ll2xy.F lnblk.F mktby_data.F mktby_parm.F my_handler.F no_digit.F op.F op1.F opencdf.F psibar.F psibavg.F psibdy.F put_bry.F put_init.F putdiag.F read1cdf.F read_grids.F read_surfnm.F readcdf.F reset_isle.F rmblklines.FF rotparm.F save.F scale_err.F set_baro.F set_land.F set_meantrc.F set_oerr.F set_pmask.F set_psi.F setbyclsd.F setbyopn.F setclsd.F setdepths.F shapmean.F spline.F splint.F state_ini.F swapdata.F transbdy.F transobdy.F trc_names.F ts2psi.F vel_err.F velocity.F vertgrid.F vertgrid_old.F wrt_transbdy.F xread.F xtr_psiby.F xy2ll.F baropar.h crs_dom.h cstseg.h curflds.h fdiag.h fne_dom.h grddat.h hybrid.h iounits.h meants.h moddat.h ndimen.h obserr.h oldflds.h param.h pconst.h pefldid.h pi_netcdf.h qgfldid.h shapfil.h switches.h version.h workspa.h zdat.h GNUmakefile.alpha GNUmakefile.cray GNUmakefile.iris GNUmakefile.rs6000 GNUmakefile.sun3 GNUmakefile.sun4 GNUmakefile.sun5 UPDATES VMSmakefile.com compile.com make_transbdy.in pe_initial.in transport.dat ============ The Makefile ============ Currently, there are seven different GNUmakefiles for seven different computer architectures: GNUmakefile.alpha GNUmakefile.cray GNUmakefile.iris GNUmakefile.rs6000 GNUmakefile.sun3 GNUmakefile.sun4 GNUmakefile.sun5 These makefiles are written for GNU Make, version 3.75 Free Software Foundation (617) 876-3296 675 Mass Ave. gnu@prep.ai.mit.edu Cambridge, MA 02139, USA http://www.gnu.ai.mit.edu/ They are NOT compatible with the standard UNIX Make. The GNU Makefiles have been designed to allow the user to compile the PE_INITIAL source codes in a separate directory from that in which the source codes are located. The makefile searches for the code segments in the following alternate paths: source code: (1) the directory containing the GNUmakefile. (2) the directory specified by the macro SRCDIR include files: (1) the directory containing the GNUmakefile. (2) the directory specified by the macro PARAMDIR (3) the directory specified by the macro SRCDIR (4) the directory specified by the macro NCDIR This provides the user with the flexibility for the following configurations: (1) The user needs only copies of this GNUmakefile, the parameter file "param.h" and a path to the source codes to produce a version of PE_INITIAL with the appropriate C-preprocessing and compiler options. (2) The user who is modifying the PE_INITIAL code, can isolate those routines actually being changed with a copy of this GNUmakefile in a sub-directory. For VAX/VMS systems, two scripts are provided for C-preprocessing and compiling: VMSmakefile.com (driver script) compile.com These are simple scripts which assume everything is in one directory. They have not been tested in years, use with caution. ------------------------------------- GNUmakefile Tunable macro definitions ------------------------------------- The User needs to check and modify the following macro definitions in the appropriate GNUmakefile before compiling and linking the application code: BIN name of the executable code. BINDIR directory path for the executable code. CPPFLAGS C-preprocessing flags and options. FFLAGS flags to the fortran compiler. NCDIR directory path for NetCDF include files PARAMDIR alternate directory path for the include files. SRCDIR directory path for the source code. ----------------------------------- GNUmakefile C-preprocessing Options ----------------------------------- The following are the available C-preprocessing options to use in the macro definition CPPFLAGS: aixdate AIX intrinsic date routine (IBM RS6000). craydate CRAY's intrinsic date routine. decdate DEC's VMS intrinsic date routine. gendbg Generic Debugging: Preserves intermediate files leqstate Linear Equation of State reset_h Reset topography for new fluxes (PE 8.3 and newer) rmcomments Remove comments when pre-processing. rmdocinc Remove documentation in all include files rmpsimean removal of level by level dynamic height mean value sundate SUN's intrinsic date routine. sunflush regularly flush output buffers in SUN systems. sunfpe enable SUN's Floating Point Exception trap. unesco Unesco (1980) nonlinear equation of state for seawater Temporary (experimental) options: gridold Use old algorithm to set vertical coordinates (backwards compatibility). meantracer Compute mean tracers for shapiro mean removal trytrby Try imposing transport streamfunction BC's in open domains. Use transport.dat. pepsib Compute PSIBAR with PE model weighting ubpsib Compute PSIBAR with weighting for uniform UBARO Temporary (experimental) biological options: bioDuse Dusenberry biological model trcascii ASCII file input for biological tracers readbioerr read error fields for biological tracers meantracer_zip Set mean tracer to initial tracer values ------------------------------- GNUmakefile Installation Issues ------------------------------- A number of internal macros are defined for the system commands used by the GNU makefiles. These will generally only have to be defined once, the first time the user installs PE_INITIAL on a new system. RMBLKLINES The name given to executable code (provided with this package) to remove blank lines from the pre-processed code. This is provided only to avoid possible conflicts. SHELL The shell to be used by the makefile. RM The remove command. ECHO The echo command. FALSE The false command. Does nothing, but does it poorly. (raises error flag) LIB The netCDF library CPP The C Pre-Processor. FC The FORTRAN compiler ===================== Compiling and Linking ===================== Once that the software has been installed and the Makefile has been selected and customized, the User needs to attend the following steps to compile and link the application code: (1) Customize include parameter file "param.h". The User needs to set the following parameters: MCSEG maximum number of coastal segments. MISLE maximum number of islands in the model basin. MPROF maximum number of points in mean TS profile. MPSEG maximum number of points per coastal segment. MT maximum number of tracers. Minimum is two for temperature and salinity. MX maximum number of points in any of the two horizontal dimensions (MX=max(IM,JM)). NP maximum number of horizontal points (NP=IM*JM). NZ maximum number of vertical levels. It is recommended to set these values big to avoid recompiling each time, for example: parameter (np=60000,mx=300,mt=2,nz=20,mcseg=20,mpseg=2000, * misle=20,mprof=5000) ================== Input Script Files ================== ------------- pe_initial.in ------------- There are several parameters that the User needs to determine and specify before running the application code. These parameters are read from standard input. A sample of the input parameter data is provided in script file "pe_initial.in". The input parameters and switches that the User needs to determine are: Card 1: IFILE Type or origin of input data: [0] NetCDF file. [1] objective analysis (OA) or feature models (FM) hexadecimal file. [2] quasi-geostrophic (QG) model hexadecimal file. Card 2: JOB Job type: [0] generate initial and boundary condition fields. [1] generate assimilation fields. [2] generate initial and boundary condition fields using the barotropic mode data provided in file FNAME(8) (see Card 27). [3] initial and boundary conditions with horizontally uniform T,S and zero velocity. Used to obtain pressure gradient and/or shapiro filter biases after removing mean T,S. *** Note: This option supersedes IFILE, the input file (3rd line file, Card 22) is a single profile (same format as 4th line file, Card 23). Card 3: TSTART Day of the input data that is assigned for initialization or assimilation. NOTE: for OA/FM hexadecimal input, TSTART and TSTOP (Card 4) is the modified Julian day (for instance, 8920) or the year day (for instance, 125) whereas for QG model input it is the model day (for instance, 1.0). The modified Julian day is defined to have an offset of 2440000 from the standard Julian day. Card 4: TSTOP Last input day to process. For persistent boundary conditions set TSTOP=TSTART. Card 5: TSKIP Time interval (days) for scanning input data. If TSKIP>0, output fields having a time-coordinate, like initial conditions, boundary conditions and assimilation fields will be processed and written every other TSKIP days. At input, the fields are needed every TSKIP days. However, if TSKIP=0, then the input data is scanned and processed at its own (input) time frequency. Card 6: ITS2PSI Flag for computing streamfunction from the TS fields.: [0] No, streamfunction or dynamic height is read from the input files. [1] Yes, the streamfunction is computed internally from the supplied T/S fields. Card 7: IBRY Type of boundary condition on transport streamfunction: [1] specify streamfunction at the boundaries and the next interior point. [2] specify transport streamfunction and the time-rate of change of barotropic vorticity (CFvN). Card 8: IFLAG(5) Switch for scheme used to compute barotropic velocity: ------------------------------------------------------------------- | reference surface | baro velocity | bdy values ------------------------------------------------------------------- 0 | take as is | do not reset [1] | QG psi [2] 1 | take as is | do not reset [1] | ubaro estimate [3] 2 | take as is | reset [4] | QG psi[2] 3 | take as is | reset [4] | ubaro estimate [3] 4 |set to DEF_DEPTH [5]| do not reset [1] | QG psi [2] 5 | set to DEF_DEPTH | do not reset [1] | ubaro estimate [3] 6 | set to DEF_DEPTH | reset [4] | QG psi[2] 7 | set to DEF_DEPTH | reset [4] | ubaro estimate [3] 8 |set from nc-file [5]| do not reset [1] | QG psi [2] 9 | set from nc-file | do not reset [1] | ubaro estimate [3] 10 | set from nc-file | reset [4] | QG psi[2] 11 | set from nc-file | reset [4] | ubaro estimate [3] ------------------------------------------------------------------- Notes: [1] set guess for barotropic velocity u_baro as estimated from psi field (u_QG) and invert for vorticity of u_baro to obtain transport. [2] bdy values obtained by vertical integration of QG psi. [3] bdy values obtained by integration of u_baro along the bdy. [4] u_baro=u_QG - a(WBOT_DEPTH) u_BOT.grad(h)/|grad(h)| (see input Card 10) if WBOT_DEPTH > 0, a=1 if h WBOT_DEPTH. if WBOT_DEPTH = 0., a=1. [5] see input Card 9 [6] see input Card 28 Card 9: DEF_DEPTH Default reference depth (m). (See description for card 8) Card 10: WBOT_DEPTH Bottom steering depth (m). (See description for card 8) Card 11: IFLAG(6) Switch for ageostrophic motion adjustment: [1] none. [2] beta effect. [3] full Coriolis. Card 12: NT Number of tracer to process. Minimum two for temperature salinity. Card 13: IBIOM Biological model type: [0] none. [1] McGillicuddy et al. (1995). Four biological tracers: phytoplankton, zooplankton, nitrate, and ammonium [2] Anderson expansion of McGillicuddy model. Five biological tracers: phytoplankton, zooplankton, nitrate, ammonium and detritus. [3] Fasham et al. (1990). Seven biological tracers: phytoplankton, zooplankton, nitrate, particulate organic nitrogen, dissolved organic nitrogen, and bacteria. Card 14: ITR Switch to activate the reading of a another input file which contains the data for additional tracers (like biological tracers): [0] no. [1] yes. Card 15: INTOPT Option to interpolate fields to terrain-following coordinates: [0] linear interpolation. [1] cubic spline interpolation. Card 16: IFILTER Switch to Shapiro filter fields: [0] filter nothing. [1] filter tracers on flat levels [2] filter velocities on flat levels [3] 1 & 2 Card 17: NORD Order of Shapiro filter (2, 4, or 8). Card 18: NREP Number of Shapiro filter applications. Card 19: IDBUG Switch to write out fields for debugging and browsing: [0] no. [1] yes. Card 20: FNAME(1) Output NetCDF file name for initial and boundary fields, or assimilation fields to be read by the PE model. Card 21: FNAME(2) Output echo (history) file name which includes diagnostics for the fields processed. Card 22: FNAME(3) Input NetCDF or hexadecimal file name to be processed. If JOB=3 (card 2), this is an ASCII file containing a single T-S profile (like card 23). Card 23: FNAME(4) Input file name for mean T-S profile. This data is used by the PE model to compute the mean density. Card 24: FNAME(5) Input GRIDS NetCDF file name. Card 25: FNAME(6) Input file name which contains information about the transport around islands and along open coastal segments. If the GRIDS file contains Land/Sea mask data, the package reads the contains of this file. Card 26: FNAME(9) Input file name which contains the additional tracer fields. The time clock on these tracer fields MUST be the same as in the primitive (temperature, salinity, velocity) fields. This file is read only when ITRC=1. Card 27: FNAME(8) Input file name which contains barotropic velocity fields and transport boundary conditions. This file is read only when JOB=2. Card 28: FNAME(9) Input NetCDF file name which contains a reference surface for dynamic height calculations. (See description of Card 8.) ---------------- make_transbdy.in ---------------- There are four files that the user needs to specify in order to extract boundary conditions for the transport streamfunction. Two file define the geometries of the respective domains, the third provides the transport streamfunction data while the fourth contains the results. Card 1: CRSGRID Name of input coarse domain GRIDS netCDF file. Card 2: FNEGRID Name of input fine domain GRIDS netCDF file. Card 3: CRSPEINI Name of input coarse domain PE_initial netCDF file. Card 4: TRNSFILE Name for output transport streamfunction boundary condition file. ==================================== Preparation of the GRIDS NetCDF File ==================================== The bottom topography for any particular domain needs to be processed before it is used by the "pe_initial" package and later by the PE model. It needs to be smoothed and conditioned to satisfy the hydrostatic condition criterion. The grid preparation package GRIDS MUST be used to extract and manipulate the bottom topography and to set up all the relevant parameters for the terrain-following vertical coordinate transformation. For more details refer to the GRIDS User's guide located in "pub/HOPS/Grids". ===================================== Preparation of the Land/Sea mask data ===================================== If the model basin includes coastlines and islands, the User needs to set-up the Land/Sea mask data which include the mask at the tracer and velocity points (recall that the PE model has a staggered B-grid), the coastline segments, and the corners of the islands boxes (if any). These mask data is stored in the GRIDS NetCDF file. The package PE_MASK appends the mask data into the GRIDS NetCDF file. This mask package is a set of FORTRAN routines and MATLAB M-files that can be used to create the mask data. Both the FORTRAN and MATLAB routines perform the same task. The MATLAB routines have the advantage that the mask data can be determined graphically. For more details refer to the file "Readme.mask" located in "pub/HOPS/PE_mask". ========= Executing ========= The input data for PE_INITIAL can have different number of levels and different depths than those of the PE model application. Once that the PE topography-following grid is determined, the initial and update PE model fields are obtained via linear or cubic-spline interpolation from the input fields. Since extrapolation is not desirable, the shallowest and the deepest input data "flat" levels must bound the shallowest and deepest PE model grid point. To execute the PE_INITIAL program, use the following command UNIX command: pe_initial < pe_initial.in > & log & To execute the MAKE_TRANSBDY program, use the following command UNIX command: make_transbdy < make_transbdy.in > & make_transbdy.log & ================= Diagnostic Output ================= One of the major jobs of PE_initial is the construction of velocity estimates from either a supplied streamfuncion/dynamic height field or from the temperature & salinity fields. To provide some quality control, PE_initial now outputs some velocity fields not directly used by the PE model. (Note: For the static initialization JOB=3, the velocity fields are trivial and these intermediate fields are not produced.) These fields are in the netCDF file under the following names: vtot - Final total velocity. vbaro - Final barotropic velocity. vgeo - Total geostrophic velocity interpolated to model levels. When the variable IDBUG is nonzero, the following intermediate velocity fields are included: vflt - Total geostrophic velocity on flat data levels. vb1st - First guess to the barotropic velocity, the vertical average of "vgeo". vb2nd - Topography adjusted guess to the barotropic velocity, only if the variable IFLAG(5) is set to reset the barotropic velocity. As an aid to interpreting these fields, an outline of the order of operations used in PE_initial to create the above fields is presented. (1) Using geostrophy, an estimate of the total velocity, VFLT, is constructed on the original flat data levels. This field should be largely horizontally non-divergent, ignoring the effect of Coriolis. (2) The results of (1) are interpolated onto the model terrain-following grid, producing VGEO. The interpolation process may add horizontal divergence. (3) The results of (2) are vertically averaged to produce the first guess to the barotropic velocity, VB1ST. (4) The results of (3) are subtracted from the results of (2), leaving the baroclinic velocity, VCLIN. This field is directly used by the PE model. (5) If so requested by IFLAG(5), the results of (3) are further modified for bottom steering. The component of total velocity at the bottom which is parallel to the topographic gradient is subtracted from (3) producing VB2ND. (6) The curl operator is applied to the results of either (3) or (5), producing the barotropic vorticity. This is inverted through a Poisson equation to construct a transport streamfunction, PBAR. This ensures that the net transport is non-divergent, as required by the rigid-lid approximation. This field is directly used by the PE model. (7) The results of (6) are differentiated to produce the final barotropic velocity, VBARO. (8) The results of (4) and (7) are added to produce the final total velocity, VTOT. ======== Examples ======== ----------------- Massachusetts Bay ----------------- The tar file "Ex_peini_9.10.tar.gz" contains an example which illustrates how the PE_INITIAL package is used. The Example is in Massachusetts Bay (MassBay) and uses data collected by Drs George B. Gardner (University of Massachusetts) and W. Rockwell Geyer (Woods Hole Oceanographic Institution) during a cruise on 24-26 July 1990, funded under the Massachusetts Bays Program. The input fields were created in the global objective analysis (OAG) package. In fact, they're the output from the example in that package. The OAG package can be found in "pub/HOPS/OAG". The topography file was created by the GRIDS package ("pub/HOPS/GRIDS"), given a land mask by the land masking package ("pub/HOPS/PE_mask") and had it's topography conditioned by the topography conditioning package ("pub/HOPS/Cond_Topo"). This topography file was taken from the output of the example in the conditioning package. This example also shows how the plotting package is used to produce horizontal contours and cross-sections. This plotting package can be found in "pub/HOPS/Plot". This example includes the following files: GNUmakefile.sun5 GNUmakefile used to create the executable. param.h PE_INITIAL main parameter statement. pe_initial PE_INITIAL SUN Ultra executable. grids_massbay.nc MassBay input GRIDS NetCDF file. oag_massbay.nc MassBay input OAG NetCDF file. ts_mean.dat Mean TS profile based on oag_massbay.nc transport_def.dat Boundary specification of transport (default) transport_p05.dat Controlling transport through boundaries. pi_ini.in input script used to create model initialization. pi_setbias.in input script used to create initialization for determination of bias in pressure gradient. pi_massbay_ini.echo output echo file from initialization. pi_massbay_ini.nc MassBay initialization netCDF file. pi_massbay_setbias.echo output echo file from initialization for determination of bias in pressure gradient. pi_massbay_setbias.nc MassBay netCDF file: initialization for determination of bias in pressure gradient. gmeta.cnt NCAR's plot file of MassBay horizontal contours. gmeta.sec NCAR's plot file of MassBay cross-sections. pi_ccnt.in input script to plotting package CCNT_NCAR. pi_csec.in input script to plotting package CSEC_NCAR. gs1.pal.white color palette used in horizontal contours. gs2.pal.white color palette used in cross-sections. plot_pe.gap.dat plotting package parameters for the MassBay. gulf_cst.dat coastlines/islands data from CIA dataset. sec_line.dat file used to indicate cross-section location on horizontal contour plots. pi_csec.log output of CCNT_NCAR, use to trace version number for plotting package. pi_ccnt.log output of CSEC_NCAR, use to trace version number for plotting package.