| Back to HOW_TO  |  Back to General Information  |



This document contains additions and changes to otemt1.

New forcing options:

In order to force the model by realistic pressure and wind fields, the following needs to be done:
  1. Obtain pressure fields (unpack files from our archive)
  2. Instruct otemw1 to use them.
General remarks:
Simultaneous excitation with tides, pressure, and wind is supposed to yield the most realistic account of time-dependent sea level and crustal loading. The major drawback in the modelling is incurred by open boundaries; there, the dynamic pressure response cannot be modelled. We can compute only the limiting cases of either perfect equilibrium response or  rigid response.  (One could conceive of a delayed pressure response, but that's futuristic music.) The model can be expected to work best in a setting of a well-enclosed basin and a boundary far away in the deep sea.

It is advisable to focus on the time-series that can be produced with otemw1, rather than harmonic solutions. For instance, the difference between the tide generating potential and the modelled elevation using pressure and wind forcing without tides is indicative of the dynamics of the sea basin (i.e. to what extent resonances could be excited by meteorological forcing). The loading effects should be compared to the bracketing cases of either pressure acting on land only (i.e. assumption of a perfect invcerse barometer response of the basin) or loading pressure everywhere ("frozen water" response).

Once this routine has computed the dynamic water loading effects you need to add the static part. Loading over land globally can be computed with the programs in  /home/hgs/Oload/smhi , aploading, or higher level shell scripts like  run_aploading_4_{project}. What's missing is the static pressure effect on the sea bottom. This is computed with another program in the OTEQ family,  apload.

To restore files from our archive, copy a monthly file

   cd /scratch/hgs/
   cp /GNSS_archive/misc/hgs2/SMHI/GE{period}.asc.gz .

and unpack
   ~hgs/Oload/smhi/aploading -p year month {month ...}

This generates files with names like /scratch/hgs/gap1997/m_oct/d31.18 from   GE9701.asc  (example is for UTC 1997-OCT-31 18:00:00). The files come every 6 hours. They are binary; subroutine GETRM can be used to input them into the processing programs. The field variable is sea level pressure minus 1000 hPa.

Currently, the file copying is through computer asmund:
   cd /scratch/hgs/
   ftp asmund userid password
   cd /GNSS_archive/misc/hgs2/SMHI/
   get GE{period}.asc.gz

Namelist &param - differences with respect to otemt1:

Obsolete parameters:
QOut_HarmQHZ, QHM, QHDC  (no second harmonic option)
PHACUT   (implementation desirable)
IDTTEP  (changed timing scheme)
QSETEL   (no start from an elevated state)
QProc_Buoy, QSafe_Buoy, QBuoy, IUBuoy  (no bouys; would be a nice feature though!)

Time intervals:
DTTGP, DTMFF - interval [s] for the recomputation of tide potential and meteorological fields,
respectively. Delfaults:  600 resp. 3600 seconds.

Pressure modelling:
QWIND, WindC, WindP(10), IWtyp, Qdo_AP, Qdo_WIND

QAP          - logical - .true.: Use the following pressure parameters:
IAPT         - integer - Model type:
                          0 - PM0 - Travelling Gaussian bell (APC,APV*,APW*)
                          1 - PM1 - Gradient East (APC)
                          2 - PM2 - Gradient North (APC)
                       > 10 = file unit number - PMF, global air
                              pressure fields  (APC,APPATH).
APC          - real    - PM0: Centre pressure in terms of inverted barometer
                              effect, given in units of meter water column.
                         PMF: Bias to be removed from file data. The effect
                              of 1000 hPa is already removed when unpacking
                              the files. The remaining bias concerns the
                              ellipsoidal component of global pressure.
                              Since it's the temporal variation of pressure
                              that drives our model, the ellipsoidal part
                              can be masked away efficiently with a constant.
                              It's on the order of 0.1 m.
APX,APY      - real    - PM0: Starting position of pressure anomaly, in grid
                              units. It's are favourably outside the area.
APVX,APVY    - real    - PMO: Velocity (km/h) east (x) and north (y).
APWX,APWY    - real    - PMO: Half-width (km).
APPATH       - char*64 - PMF: Symbolic path to the unpacked pressure files.
                              Default is APPATH=' ', i.e. the built-in default
                              '/scratch/hgs/gap%%%%/m_@@@/d##.## '
                              Replacement: '%%%%' by year, '@@' by decimal month
                              '@@@' by lower-case character month,
                              '$$$' by upper-case character month, first
                              occurrence of '##' by day and second '##' by hour.

QAB_EAPR     - logical - Equilibrium pressure response enforced on the open
Qdo_AP       - logical - Pursue pressure modelling.
QAPSTOP      - logical - End pressure forcing .
OPT_QUINT    - char*16 - PMF: Send an option to interpolation routine.
                              This could in the first case be the trace
                              option, '+T'.

QWIND        - logical - Use the following wind parameters:
IWtype       - integer - Model type
                         0,1,2 - WM0 WM1 WM2 - corresponds pressure model  (IAPT)
                          >= 3 - WM3  - compute geosptrophic wind from pressure.
             mod(IWype,10) = 4 - WM3' - set pressure field to zero after
                                        having computed the wind field).
WindC        - real    - PM0: Centre pressure. Gaussian bell anomaly.
                              Wind is derived accordingly.
                              WindC should be the same centre pressure as APC.
                         PM1, PM2: Wind shear. The supplied value is to represent
                              sigma*DTN/RHOW (times time-step devided by density
                              of sea water.
WindP(10)    - real    - WM3: Geostrophic wind parameters. See below.
Qdo_WIND     - logical - Pursue wind modelling.

OPT_AP_EVWH  - char*1  - PM0 PM1 PM2:
                              'Y' - Fill the area with pressure effect
                                    at all nodes
                              'N' - ...only at sea and active-boundary nodes.
                              ' ' - internal default = 'N'
                         WM3: ' ' - internal default = 'Y'

Wind model and geostrophic wind parameters

            wind'=(WindP(1)+WindP(2)/|grad p|)*[grad p]/(f rho)

            wind=Rot(WindP(3))*wind'             rotation angle in degrees
            stress=WindP(5)*wind*|wind|          if |wind| < WindP(4)
            stress=WindP(9)*wind*|wind|          if |wind| > WindP(8)

            Rot(a) is a rotation matrix for angle a
            [v] is a 90 degree rotation of vector v
            |v| is absolute value of vector v

            Default values are used unless supplied value =/= -999.9

                  WindP/ 0.56,    2.4,        20.0
                          5.0,    0.565e-3
                     -0.12e-3,    0.137e-3
                         19.2,    2.513e-3,    0.0/

To execute, issue
   cd ~/OTEQ
   otemw1 MODEL_SUBDIR/instruction_file[.ins]

 See an example instruction file: MDL22h/otemwww.ins