============================ OCEAN LOADING POST PROCESSOR ============================ INSTRUCTIONS FOR USE Hans-Georg Scherneck Chalmers University of Technology Onsala Space Observatory July 2, 1997 1. Background A short account on the loading problem: The load-induced displacement due to an oscillating point load with mass dm is du = Diag[1, sin(azi(r_load,r_site)), cos(azi(r_load,r_site))] * * [Green(radial; angle(r_load,r_site), Green(tangential; angle(r_load,r_site), Green(tangential; angle(r_load,r_site)] * dm where r_site = the position of the site at the earth surface..., azi = azimuth and angle = distance angle, assuming a sphere (ellipsoidal shape of the earth can be neglected). du = 3-D displacement = [u_vertical, u_west, u_south] = Real(U exp(i 2pi f)) dm = rho z ds z = Real(Z exp(i 2pi f)) Z (complex) specifies the amplitude and the phase of a partial tide (frequency f), rho = water density, and ds is a surface element. The tide phase is taken with respect to the luni-solar tide at zero meridian; positive phase values specify a lag (oceanographer's convention). U (complex) is the interesting quantity (3-D displacement amplitudes and phases). Island sites in the Pacific might experience up to 10 cm radial displacement. Tangential displacements are usually 1/3 of the radial ones, with a tendency to disappear on small islands and to increase towards the coast of large continents. In most cases the total displacement is on the order of 1...3 cm. Global tide models specify in tabulated form Z on a global grid (Schwiderski's geographic resolution = 1x1 degree). Programs OLM613 or OLFG do the global computation: They step through a table of a partial tide and convolve Z(longitude,latitude) with the Green's functions for the radial and tangential displacement components. A set of sites can be processed simultaneously provided the particular Green's function applies. Thus, a set of inland sites can be processed together or a set of island sites, respectively. The global load step needs to be followed up by a local refinement step if the site is near the coast (100 km); this is the task for OLMPP. Don't forget that the programs will only compute the loading effects due to the (leading) eleven partial tides covered by the ocean model. The next five hundred partial tides contribute roughly one tenth of the effect. If linear interpolation of admittance coefficients (loading effect / tide generating potential) in the frequency domain were permissible, you could arrive at millimeter precision. However, I am hesitant since there is nonlinearity in the shallow-water response. OLM613 has been run successfully on an IBM compatible PC under MicroWay's NDP-Fortran77 and on an IBM mainframe under MVS during the '80s. OLFG, originally running on the same PC, has been adapted for UNIX and has been extended to prepare exception tables - the data passed to OLMPP. In the following, the first computation stage of OLM613 or OLFG are equivalent; OLM613 is synonymous for OLFG. However, OLM613 has never been brought to UNIX. However, there is one important difference, an obstacle which is to be knocked off in the next revision: OLFG does not retrieve TOPO data. Instead this is deferred to the OLMPP stage. The solutions have been checked in several ways (comparison with Goad & Agnew's computed values for the Mojave site; comparison with observations for the case of gravity perturbations; Schuh & Moehlmann (1989) found some reduction of residual error in VLBI). NGS, IERS and GSFC use my figures. As regards ocean loading effects on SLR and VLBI the matter seems a little premature (so don't be overoptimistic as regards GPS). Henrico Derks, then at Dept. Aerospace Engineering at Delft University, NL, has done some analysis with SLR. Little of the effect could be seen although he solved for tidal Love numbers. I have his report here; you may ask Henrico, his former department (Ron Noomen), or me for a copy. During the next years we will see whether observation is in support of my computations. As regards gravimetry the results are excellent (thanks mostly to Schwiderski; his models appear to work best in application to loading problems); that conclusion is based on a comparison with O. Francis' results (he is at the Royal Observatory, Bruxelles) The operation of OLM613 (OLFG) is of type "batch-execution". Instructions for use are contained in the commentary of file OLM613.F (olfg.doc) OLMPP works interactively, using graphics. The bulk of this manual will help you through the various phases of OLMPP. The program computes the loading effects from ocean tides inside an area of 3 x 3 deg around a site. It post-processes the results of OLM613 on the basis of the "Exception Table" produced with the X-flag set. The near-field effects computed by OLMPP are added to the far-field effects of OLM613. The resulting, completed loading coefficients are written to a file. In the revision of autumn 1994, OLFG prepares data in a 3 x 3 degree area cut out from an ocean tide model while handling grid parameters more flexibly. Now, also Le Provost's 0.5 x 0.5 deg model can be processed. Do you really need OLMPP ? If you need an accuracy of 1 mm (for the total effect of the produced ocean loading time series), and if your site is near the coast, and if the tides are high there you might need OLMPP. If your site is on an island so small that it is not represented as land in the tide map, you should use OLMPP. However, if some-one claims that Schwiderski's tides are too inaccurate and the partial tide approach poses a fundamental limitation at this level, i.e. that an accuracy of better than 5 mm is doubtful anyway, I would not strongly oppose. My only counter-argument would refer to experience with tidal gravimetry; there, the closure error (solid earth tides + ocean loading - observed effect) is generally 5 percent of the loading effect. A radial displacement of 1 cm would translate into an error of 0.5 mm for a single partial tide. The "improvement" of the modelled loading effect due to OLMPP is typically on the millimeter level. Since my experience is more or less limited to Western Europe, I might be proven wrong in making this general statement. Post-processing does not have to be considered if the site is 100 km or more inland. If the land-sea resolution / representation in Schwiderski's maps is very bad near a site, your results will not become better with the post-processor. You should try to obtain a regional tide model instead. Contact me if you need a programme for the loading effects of a regional model. Such models can appear in so many different kinds of formats and grids that a general processing facility appears just impossible to code. The Exception Table output by OLM613 contains for each site: Global loading results Near field tides (TIDE) Near field ocean/land flags at a resolution of 12 per degree (TOPO) OLMPP inputs the Exception Table. A high-resolution flag array, called OCEAN, is used to refine the ocean/land distribution interactively. The initial version of OCEAN is contructed from TOPO. The geographic resolution of OCEAN can be set to 1 ... 4 times of that of TOPO. TOPO and OCEAN cover an area of 3 x 3 degrees, TIDE an area of 7 x 7 degrees. Can you use OLMPP at all ? There is a version for 386 or 486 PC's with DOS-extenders, requiring the special compiler of late MicroWay. You also need a VGA video card that supports screen mode 18 (640x480 pixels, 16 colors). You will need the LIBGREX.LIB library of MicroWay, their Fortran-77 compiler, and linkage tools to assemble protected mode executables (e.g. PharLap). On an AT 386 Personal system, OLMPP can be executed under e.g. PharLap's DOS-extender RUN386.exe Program module OLMPP.exp contains protected mode code. It requires the 387 coprocessor and 1 Mbyte of memory. Make sure the mathematical subroutine library contains REAL*8 FUNCTION DARSH (X), the inverse hyperbolical sine of real*8 X. There is another version for UNIX using the pgplot library for screen graphics under Xwindows. You need to update the pgplot library (nasty) and install a set of interface routines (4mw) to mimic MicroWay's libgrex. That version runs on HP-UX Version 9. No other experience. 1.1 Computation of the near-field loading effects. The tides are linearly interpolated (c.f. Chap. 1.2) on the OCEAN grid from the values available in the 7x7 TIDE grid. The loading kernel function is computed for each of the elements of OCEAN. Kernel function = G(distance)/(2 sin distance/2) For computational speed G is tabulated in a long array. The range of the functional values of G is limited (between 5 and -1); the variations of G w.r.t. distance are evenly distributed if distance is mapped on a loga- rithmic scale. The user has limited influence on these features. He/she, however, has to Select a Green's function file according to the structure of crust and mantle in the area of the site (oceanic, standard continental, thick shield). The Green's function currently in effect will apply to all sites processed in one OLMPP session. The user can invoke another Green's function only by specifying it in the Namelist-block &PARAM. This should not imply a limitation, since the same convention holds also for each execution of OLM613 and OLFG. 1.1.1 Mass conservation at the coast. The reason why the mass conservation feature was introduced into OLMPP is that Schwiderski NOT ALWAYS specifies observed tide AMPLITUDES whereas the loading effects are ALWAYS due to tide MASS. The central question is: What would happen to the hydrodynamic tide solution if the coastline resolution increases ? Imposing mass conservation everywhere by employing the equation of continuity, tide amplitudes tend to increase as land grows into a box of the 1 x 1 degree model, and decrease as land retreats. However, this is not a strictly locally limited phenomenon. If the area in question s excited near a normal mode, the amplitudes in a wide area may change (best example: Canadian tide energy plant building consequences in the Bay of Fundy). The 3 x 3 degree area does not strictly conserve its mass, of course. Consequence for loading effects: The impact of mass conservation is maximum if the masses are redistributed in the vicinity of the changed coastline. OLMPP will redistribute the masses inside specifically flagged 1 x 1 boxes. Since there are no hydrodynamics involved in that computation, the mass is regarded as a uniform slab. This is a middle-road between a scaling-up of the land-blanked mass inside the same TIDE box, and a distribution of the mass into an area exceeding the limits of the 3 x 3 (7 x 7) degree region. On the other hand, if the tide model would specify observed amplitudes every- where, we simply have to blank out those OCEAN elements which contain land. llllllllllllllllllllllllllllll llllllllllllllllllllllllllllll llllllllllllllllllllllllllllll llllllllllllllllllllllllllllll llllllllllllllllllllllllllllll __ llllllllllllllllllllllllllllll __ llllll..................llllll llllll......LL......LLLLllllll llllll..................llllll llllll................LLllllll llllll..................llllll llllllLL................llllll llllll..................llllll llllllL.................llllll llllll..................llllll llllll.................Lllllll llllll..................llllll llllll.................Lllllll llllll..................llllll __ lllll,..................,lllll __ llllll........................ lllll,........................ llllll........................ llll,,........................ llllll........................ llll,,........................ | | | | 1x1 deg resolution High-resolution L - "new" land-elements , - "new" ocean-elements If tide AMPLITUDE is known here, loading masses are removed as land elements are introduced. If tide MASS (VOLUME) is subject to the equation of continuity, the L-elements tilt the mass balance. The mass has to go somewhere. The ,-elements tilt the mass balance in reverse: the mass has to come from somewhere. Fig. 1: The two possible consequences of introducing Land and Ocean elements in the high-resolution OCEAN grid. In the cases where ocean-covered elements are located in a land-box of TIDE, a tide value must be inferred. OLMPP does that by first interpolating between neighbour boxes that share an edge; if this fails an attempt is made to extrapolate from all the other TIDE boxes. The inferred value can be quite weak. The user is NOT encouraged to introduce sea covered OCEAN elements if the nearest available tide box is remote. Schwiderski's model is a mixture of the above extremes: It specifies observed amplitudes at places where tide gauge information has been adopted. Inbetween, the tides are hydrodynamically interpolated; thus, some coastal boxes conserve mass, but most of them don't. Le Provost's or Ulrike Seiler's models are not constrained by prescribed tide height. We might interpret these tides as WATER MASS TIDES everywhere. It will be the aim to initiate default mass conservation flags that best fit the kind of tide model being processed. At the OLM613 stage we don't make any distinction between the two cases - their effect is insignificant (I've checked that by extending the above procedure to the global case, assuming the extreme case: no mass conservation anywhere, and computing loading effects for the most sensitive instrument: gravimeters). Near an observing site, distinction between the two cases, which are called Constrained - Tide gauge data has fixed the amplitude Unconstrained - Mass conservation may alter the amplitude after OCEAN update might be important. OLMPP handles them using either strategy described above. In order to check the impact of coastline refinement on the computed loading coefficients, the user obtains - roughly speaking - the two ends of the error bar by: (1) Recomputing the loading effect in "Single-point" mode while refining the coastlines; output the result. (2) Assigning Mass Redistribution flags (MRD(i,j)=.true.) and recomputing the loading effect of the entire area (g-key); output the result. After reading Chapter 2, the missing pieces will have fallen into place. 1.2 Interpolation Tide values from Schwiderski's maps are on a 1x1 degree grid (LRG). We are confronted with the problem to determine tide values on our high- resolution grid (HRG). There are a number of possibilities. The most advanced one would be to solve the hydrodynamic equations on the HRG, where excitation would be on remote open boundaries. This appears too difficult. The simplest assumption is: constant values inside a LRG cell. Interpolation is difficult since LRG values on land are undefined; i.e. we have a sparse suspension array. There are two interpolation routines included in the program package. OLINTP.F does linear interpolation. The three triangles of suspension points with vertices closest to the point of interpolation are detected; gradients of tide values in the three triangles are computed and their mean is used to interpolate. The suspension array consists of a 5x5 LRG. It is copied from the 7x7 tide grid. As the number of land points on the LRG increases, the number of triangles might become less than three. In extreme cases, no triangle can be found in which case a constant tide value must be resorted to. OLINTPG.F does polynomial interpolation; the degree of the polynomial is decided upon using a generalized inverse approach. A priori weights must be specified for the nearest tide value and a distance dependent weight function parameter for the nearest-but-one neighbours. The suspension grid is the same as above. The generalized inverse yields a parabola even in the case that only two suspension points are available. In order to control the curvature of the parabola, a "hydrodynamic constraint" was derived from the tide equations, which relates the Laplacian of the tide elevation to the frequency. Coriolis effects and other complications are neglected. Therefore, the weight with which the hydrodynamic constraint is obeyed should be low. The effect of the hydrodynamic constraint is to yield lower curvature (longer tide wavelength) for the low-frequency tides. Compiling with OLINTP: Two calls to GI5x5* should be deleted in OLMPP.F Associated parameters may be deleted from the Namelist PARAM. Compiling with OLINTPG: OLMPP calls parameter setting subroutines GI5x5S and GI5x5H (defining the a priori weights). These are input into the program via Namelist PARAM (c.f. Chap. 2.1.1) Both interpolation routines have the same name: ZINTPV. Create two library files and include the desired one in the link process. 2. Operating OLMPP 2.1 Prerequisits Before the user starts OLMPP, he/she should have (1) an "instruction file" available, named OLMPP.INS containing Namelist data. Example: File OLMPP.INS Namelist parameters' meaning c.f. 2.1.1 It is recommended to create a subdirectory for the data files. Each subdirectory would contain one working copy of OLMPP.INS. (2) Green's functions tables available under directory path and file name as indicated in the Namelist variable GREFCT. The same for the Exception Table. (3) ETOPO5 or TerrainBase data (from NOAA/NGDC) and a compatible reading subroutine. Reformatting is required, c.f. chapter 2.1.2 In OLMPP.f there is exactly one CALL RTOPO. The subroutine is in otes42.f. See chapter 2.1.2. (4) optionally, coastline polygone data. Serves rather for coarse orientation than for precision land-sea editing, but can be quite useful. The path name default is something/GLCOAST (set in DATA, variable through NAMELIST). See chapter 2.1.3. (5) geographic maps of the site areas available (scale e.g. 1:2,000,000). (6) when processing Schwiderski tides, (one of) his publications (Subtitle "Atlas of Tidal Maps and Charts"). After successful execution, OLMPP will have (1) created or appended a result file (FOUT='name') for loading coefficients; (2) stored away the screen palette on OLMPP.CLR for later restore. (3) for each site a TOPO file in the 3x3 degree area can be saved ("SITE- TOPO", "TOPO-SAFE"), a user-edited working copy. 2.1.1 Namelist &PARAM parameters The Namelist is input from the .ins-file once and only once per run of OLMPP. Parameter Type Meaning [Default values] USERNM - char*65 - A message like 'Anders Celsius, Uppsala University, 1726' to denote the origin of the coefficients in the #BLOKQ output record. WKDIR - char*64 - Working directory. Specify '?' to obtain the directory of the actual *.ins file. ['./'] QA - logical - Automatic processing. Prompts are largely supressed. Interactive update of neither TOPO nor OCEAN arrays. [.FALSE.] IUINS - integer - Only together with QA=.TRUE.: File unit number for reading the numbers of sites to be processed. 5 - screen 4 - *.INS file [5] IXIN - integer - File unit number for input exception table [31] FXIN - char*32 - its file name. Code '?' to get a prompt ['?'] GREFCT - charr*32 - file name for the Green's table. Right delimited by ']'. ['/users/hgs/Oload/Greenstb/mc00egbc.dat]'] Abbreviate the path root as '~/' to obtain the default of greensf.f which is '/users/hgs/Oload/Greenstb/' The associated file unit is 21. IOUT - integer - File unit number for resulting loading coefficients [41] FOUT - char*32 - its file name. ['olmpp.blq'] FRIN - char*32 - File name for input #BLOKQ-table, ['?' = prompt] IRIN - integer - its file unit number. [42] FFOUT - char*32 - File name for result, OLF-format ['olmpp.rsl'] IFOUT - integer - its file unit number [-1] COMMNT - char*32 - Comment on OLF-records ['_PP'] NH - integer - {1,2,3,4} High-resolution of the map, 12*NH per degree. The 12 is due to the resolution of the TOPO05 file. Every topo box is subdivided into NH x NH sub-boxes as the unit loading element. [0 = ask] PHACUT - real - Angle where to cut circle for BLQ output phases E.g. PHACUT=180.: phases -180...+180. [0.0] BLQAMF - real - Amplitude factor for BLQ output, suggest 1.0 for displacements and 10.0 for gravity. [1.0] SHORE - real - Shore level increment (meter) [10.0] LANDOFFS -integer- Default land/sea discrimination offset [-1] Tide interpolation parameters: MRDV - integer - Default value for filling the MRD flag array, indicating whether mass redistribution is allowed (1) or not (0). In conjunction with mass- conserving ocean tide models, a value of 1 is appropriate. The following default is good for Schwiderski: [0] WGHTC - real - Weight**2 for fit at center (neighbour points get a weight**2 of 1 each) [20.0] WHGTH - real - Weight**2 for hydrodynamic constraint [0.001] NEIBOR - integer - Apply hydrodynamic constraint in neighbourhood range -1: not at all; 0: only at centre; 1: neighbours) [0] ALWHYD - logical - False: Apply constraint only if less than 5 tide values are available; True: Apply constraint always. [.false.] Other: GLCPTH - char*32 - Path to coastline files. Don't code trailing "/". ' ' - keep internal default (pgm glcors.f) ']' - disable [' '] TOPO data retrieval: IUTOPO - integer - File unit number for data [71] NUTOPO - integer - IUTOPO+NUTOPO = unit number for descriptor file, -999 if all in one chunk [1] DSCFN - char*32 - Name for descriptor file (or data file name if all in one chunk) ['/users/hgs/TOPO/topo.dsc'] IUTOPS - integer - File unit number for storage of topo data [73] TOPSFP - char*32 - Path name for topo storage file ['. '] TOPSDP - char*32 - Path name for the TOPS.DIR file ['. '] MAXALIAS- integer - Max.number of searches in TOPSDP [5] If IUINS=4 and QA=.TRUE. the program expects a set of lines immediately following the namelist block. Each line contains the running number of a site selected from the exception table. The last site number record may contain a 'T' to switch back to interactive (terminal) input. QALLA - logical - Under QA=.true. (automatic site input): skip the topo update. 2.1.2 ETOPO5 data set organisation. (Substitute TerrainBase for ETOPO5) NOAA's original data consists of binary short integer (2 bytes) records without record length descriptor words (i.e. it is readable for C- or Pascal programs), whereas Fortran supposes the following record structure Record = L B1 B2 ... BL L File = Record1 Record2 ... where L is a 4 byte integer specifying the number of bytes (B) in the record. A little C-programming (Basic, Pascal, ...) will do the task. There is a little Basic program in ~/Oload/p/b/rtopo.bas The user must also find out whether the bytes must be swapped to comply to the machine's specification (HP: yes). As a trade off between speed and space, ETOPO5 has been copied in chunks of latitude rings, each 10 degrees wide. A descriptor file is used to assist the file access: -------------------------------------------------------- $ cat /users/hgs/TOPO/topo.dsc 18 0.0 4 90.0 80.0 /geo/hgs/TerrBase/Data/tbas000.dat 80.0 70.0 /geo/hgs/TerrBase/Data/tbas010.dat 70.0 60.0 /geo/hgs/TerrBase/Data/tbas020.dat 60.0 50.0 /geo/hgs/TerrBase/Data/tbas030.dat 50.0 40.0 /geo/hgs/TerrBase/Data/tbas040.dat 40.0 30.0 /geo/hgs/TerrBase/Data/tbas050.dat 30.0 20.0 /geo/hgs/TerrBase/Data/tbas060.dat 20.0 10.0 /geo/hgs/TerrBase/Data/tbas070.dat 10.0 0.0 /geo/hgs/TerrBase/Data/tbas080.dat 0.0 -10.0 /geo/hgs/TerrBase/Data/tbas090.dat -10.0 -20.0 /geo/hgs/TerrBase/Data/tbas100.dat -20.0 -30.0 /geo/hgs/TerrBase/Data/tbas110.dat -30.0 -40.0 /geo/hgs/TerrBase/Data/tbas120.dat -40.0 -50.0 /geo/hgs/TerrBase/Data/tbas130.dat -50.0 -60.0 /geo/hgs/TerrBase/Data/tbas140.dat -60.0 -70.0 /geo/hgs/TerrBase/Data/tbas150.dat -70.0 -80.0 /geo/hgs/TerrBase/Data/tbas160.dat -80.0 -90.0 /geo/hgs/TerrBase/Data/tbas170.dat ---------------------------------------------------------------- from to latitude file ---------------------------------------------------------------- Each file .../tbas###.dat is binary integer*2. One record contains a quarter of an elementary latitude ring, starting at zero and continuing eastward. This is the relevant part of otes41.f, subroutine RTOPO: integer*2 it(4320) ... do 222 irp=1,4 ie=irp*1080 ib=ie-1079 222 read (iun) (it(i),i=ib,ie) 2.1.3 Coastline polygone data This is pickpocket haul from late Disspla. Directory C:\GLCOAST shows NAMERICA BCD 183268 92-11-17 11.15 SAMERICA BCD 42764 92-11-17 11.15 AUSTRAL BCD 58028 92-11-17 13.32 AFRICA BCD 37356 92-11-17 8.40 ANTARCT BCD 16068 92-11-17 13.25 EUROPE BCD 99396 92-11-17 8.36 ASIA BCD 109384 92-11-17 8.39 GLCOAST DIR 519 92-11-17 13.32 The data files contain binary real*4. The reading routine is in GLCORFS.f 2.2 OLMPP loads data. Start processing: RUN386 OLMPP The first prompt concerns the INStruction file. You may enter (alt 1) the full path and file name, e.g. oload/gps94/olmpp1.ins (alt 2) a (sub-)directory followed by a backslash, e.g. GPS94/ The system will complete the string by adding olmpp.ins (alt 3) a Carriage Return. This gets you the default olmpp.ins in the current directory. (alt 4) You specify on the command line according to (alt 1..3) and enter a $ sign at the INStruction file prompter. The trailing .ins can always be skipped. The entries that follow are fairly obvious. There is an optional section of editing the coefficients file that was issued together with the exception table. This section is entered if FRIN is specified in the namelist. Follows the site prompt: Processing exception table file: "FXIN" Exception tables are available for the following sites: #. SITEN #. SITEN ... Specify site number to process [, [No]Auto]: The user enters the number of the site to post-process. An " F" can be appended after the site number if a new regional TOPO-grid (reFresh) is to be enforced. The next prompt asks for the desired OCEAN resolution. The number to be entered signifies the number of loading elements per TOPO element, 1..4 are legal values. After that OLMPP displays (1) the site coordinates and area location; (2) the TIDE availability flags, a 7 x 7 array containing T where values of the 1 x 1 degree tide model are available, and F where they are not (so-called Schwiderskiland); (3) the preliminary loading coefficients from the global step (OLM613 or OLFG). 2.2.1 OLMPP asks whether it shall display the TOPO-grid. Internally in OLMPP, topography is characterised by 15 height levels. The levels are -7 ... +7. The level increment in meters is entered to OLM613 through the array SHORE, SHORE(i) for site # i. For reading fresh TOPO in OLMPP, the level increment is specified through variable SHORE and the land-sea discrimination level through variable LANDOFFS. The purpose of this feature is to provide a realistic first order discrimi- nation of land and sea boxes in the local area. The default coastline level is 0. OLMPP displays the TOPO grid on a graphic screen. Topo levels below the coastline level are displayed in blue colours, above it in red colours, respectively. OLMPP prompts for a new coastline level. Enter 99 to leave the topo display, using the currently defined level. Entering positive values will raise the land-sea discrimination level (low "land" turns to "sea") and vice versa. A level is to be found which requires minimum hand work at the stage of Ocean Update (Chapt. 2.4.3), i.e. one which yields a land-sea distribution closest to reality. In flat areas, the coastline can change significantly after adjusting the discrimination level by a small step (e.g. equivalent +- 10 m in the case of Florida). Entering +8 will turn the whole area to "land", -8 to "sea". If you are unsure, don't change the default coastline level (enter 99), go on until OLMPP is in the Ocean Update stage, and judge how many of the high-resolution boxes you will have to change. If there are many boxes, judge whether there is systematic dependence on the topography. If yes, quit Ocean Update, reprocess the site, and adjust the coastline level at the TOPO-display. NEW: Individual TOPO elements can be updated. Select Update (press U-key) to enter this mode. Then press... arrow keys - to locate cursor ESC - to return to the TOPO-level prompt + - to increase topo value > - " " " " - - " decrease " " < - " " " " B - to start block mode = delineate a reactagle subject to < + > - N - to calibrate cursor position RET - to end block mode 2.3 OLMPP calls the Ocean Loading EFFect subroutine OLEFF. Next, the system calls subroutine OLEFF, which sets up / completes the arrays in its first part, and - if the user has chosen so - computes the loading effects on the current state of the OCEAN grid. Since this is the most time-consuming step, an option to skip computation on the raw grid was introduced; if a substantial part of the OCEAN must be revised, the raw grid scan is just a waste of time. At stage one of OLEFF the screen output shows the TIDE availability flags and the MRD array (Mass ReDistribution flags) in the 3x3 elements grid. Some previous False flags have turned True as OLEFF was able to infer tide information from the neighbourhood. The MRD flags are False by default; only the user can change them. If load effect computation is selected, the system refreshes a screen message at every new row of the OCEAN scan. After completion of OLEFF, the graphic section of the program is activated. The associated program section of MAIN has been broken out to form a subroutine Screen_Routines. CALL Screen_Routines can be deleted in MAIN; the result is a program that could (after a few other code changes) run as a batch job, adopting the TOPO determined ocean/land distribution ("raw OCEAN grid"). You must also delete the call to PTOPO then. The major purpose of Screen_Routines is to activate Ocean_Update for interactively adjusting land/sea distribution towards realism. 2.4 Graphic procedures - displays and user entries In Screen_Routines the system changes display mode to graphics and draws the OCEAN grid. The first routine called is Display_Results. After a displacement component has been selected for display (r - radial, w - west, s - south), control is passed to Ocean_Update. Ocean_Update supervises the cursor-prompting procedures and computational procedures. 2.4.1 Keyboard entries while screen is in graphic mode Graphic screen <----> Keyboard entries are by poll (they are not followed by Carriage Return) Keys assigned to operations: Arrow keys - 4 cursor control keys F1 F2 F3 - 3 Function keys CR - Carriage Return SP - Space key (SpaceBar) a b c ... - Upper-case characters are synonymous to lower-case < > + - ? - A few special symbols ^P - Control-P - press before you PrintScreen. The Escape-key has not been assigned; the procedures take no precautions as to be able to back up or discard recent entries. The routine Cursor_Prompt: Function keys F1 F2 and F3 are interpreted inside the procedure Cursor_Prompt. The lower right-hand corner of the screen shows that scanning of the F-keys is active; simultaneously, the cursor is located on the OCEAN display. Upon activation of a Function key, Cursor_Prompt branches internally to F1 increase, F2 decrease the cursor movement interval, F3 change the color palette. Although unnotified, ^P is checked by Cursor_Prompt to enable PrintScreen. The cursor control keys move the cursor over the OCEAN grid; other keyboard entries cause the system to leave Cursor_Prompt; Ocean_Update will interpret the entry and possibly pass control to another procedure. 2.4.1.1 Changing the cursor increment There is rare use for this option. Appropriate cursor steps and location inside picture elements are set by the program. However, the lines drawn at integer latitudes and longitudes are drawn in a color also used for the case "problem at the observing site, site box state = ocean". In order to get the cursor to odd places and temporarily change their color, cursor increment = 1 must be used (pressing F1 repeatedly). 2.4.1.2 Changing colors - temporary and permanent Concerns the palette assignement to the color code of the pixel where the cursor is located at. Entering F3 implies an internal branch in Cursor_Prompt. After F3, the cursor is hidden. The user is prompted to enter SP - to assign the next higher or lower palette color, depending on... > - to assign the next higher palette color; < - to assign the next lower palette color; - - to set the palette scan interval to 1 + - to set the palette scan interval to 8 (there are 64 different palette colors on a VGA mode 18 screen; an increment of 8 will take you through related colours); CR - to leave the branch and return to the cursor, keeping the present colors. This function is mainly for the user to distinguish between similarly colored OCEAN elements. E.g.: Discern the element which contains the observing site (the eye will find it when the elements around it change color). Every time OCEAN is redrawn, the permanent palette is taken. A permanent change of the palette setup is initiated with the p-key. The cursor is now restricted to scan the color tags of the "legend", which appears on the right-hand side of the screen. The procedure to assign palette colors is the same as above. At quit-time, the user can specify whether the new color set is saved for subsequent runs. 2.4.2 Graphic displays 2.4.2.1 Results display - one of the three displacement components This display is issued by subroutine Display_Results. The subroutine has a prompting and a no-prompting mode. The no-prompting mode is invoked during Ocean_Update editing of single OCEAN elements. Displayed are the results for the eleven partial tides (the ocean loading coefficients: amplitudes in millimetres and phases in degrees) for one of the three displacement components. The prompting mode is invoked by pressing the r-key in Ocean_Update or after OLEFF compution. Successive entries: r - selects the radial component w - selects the tangential component along west azimuth s - selects the tangential component along south azimuth CR or SP - orders to leave Display_Results. Any other key - display the currently selected component. The default component is 'r'. 2.4.2.2 OCEAN grid - states, attributes and their graphic representation The screen position where the ocean grid is drawn is fairly obvious. After the first time Draw_Ocean is called, there are the following possible states and attributes for ocean elements (pressing the x-key will order a legend on the right hand side of the screen): State Attribute Explanation land - unmarked No loading effects computed at land elements sea - computed Loads have been accounted in subroutine OLEFF sea - unmarked OLEFF was asked to skip computation site/land - unmarked The observing site is on land by default. There is only one such element. The following cases may also appear: land - permanent No possibility to infer tide information here. The user cannot edit these elements. sea - problem A TOPO element far inland has determined some OCEAN elements' state "sea". Inference of tide information was too difficult, however. The remaining cases (might) appear when editing OCEAN. Their meaning is obvi- ous. For additional information, c.f. Chapter 2.4.3, Updating the OCEAN grid. 2.4.2.3 Top right part of screen Reserved for displays of loading coefficients for one of the three components (r-key) the color legend (x-key) 2.4.2.4 Middle right part of screen The 3x3 MRD array is displayed here. Mass ReDistribution flags are valid for each 1x1 degree box of TIDE. The sections of OCEAN delineated by thin lines are represented by one flag in the corresponding place of MRD. F - FALSE = Constrained by tide gauge data = No mass redistibution here. T - TRUE = Unconstrained = Mass continuity = Mass redistribution needed. The mode of Ocean_Update is displayed here: Single-point mode (SP), Multi-point mode (MP), and whether Multi-point mode is pending. C.f. Chapter 2.4.3, Updating the OCEAN grid. 2.4.2.5 Bottom right part of screen Reserved for activated keys display. Active keys appear colored. There are too many keys allocated to edit functions to fit in this place; a help-screen is obtained after the h-key or the ?-key has been pressed. 2.4.3 Updating the OCEAN grid Help function: h or ? The h-key or ?-key under cursor_prompt requests a page of help. Help concerns the Update_Ocean edit functions and conveniences. You should see a cross-hair cursor about mid-screen. Move the cursor to an element of the OCEAN grid the state of which is to be changed; use the arrow keys. The action of the system to a state change depends on the selected operation mode. 2.4.3.1 Single-point mode (SP) Edit functions: + and - Result display: r Quit: q Editing is done on single elements of OCEAN; the change of the loading effect is computed everytime an element is changed. The cursor location determines which OCEAN element is concerned. Keyboard entries: + changes an OCEAN element from state "land" to "ocean". The additional loading effect due to this single element is computed and added to the loading coefficients. - changes an OCEAN element from state "ocean" to "land". The following will/can happen: The result display is refreshed. OCEAN is redrawn at the cursor location. The success of the computation is reflected in the color attribute of the OCEAN element. The legend helps to interpret the attribute (press x-key). The following states and attributes can appear: State Attribute Explanation Land - computed The loading effect of the element earlier in state "sea" has been subtracted. Sea - computed (reverse of the above) Sea - problem Tide information could not be obtained here. No loading effect computed. Site/sea - computed The site is now in a "sea" element. The loading effect is computed using an approximative, semi-analytic function. Site/land - computed (reverse of the above) If '-' is called for a "land" element, it remains "land". If a "land" element has had attribute "computed" it will be turned into "land - unmarked". The equivalent holds for two successive calls of "+". If tide information cannot be inferred at an element the state of which is changed from "land" to "sea", it will be marked "permanent land". r invokes Display_Results in prompting mode; the user can select another component. q quit, leave the graphic part, go to output section ... An additional key strike is necessary: s to leave Ocean_Update, saving the permanent palette; r to resume Ocean_Update; any other key to leave Ocean_Update. 2.4.3.2 Multi-point mode (MP) Edit functions: o l u c m Computation function: g l and o The user marks a number of OCEAN elements "land" (l-key) or "sea=ocean" (o-key). The elements which are going to be marked, plus the elements changed in single point mode will be considered in subsequent calls of OLEFF. During multi-point mode, OLEFF is not called at each edited element. The displayed table of load eff. is no longer valid for the present land/ocean distribution. g Computation of loading effects is ordered by pressing the g-key (Go!). The entire OCEAN grid will be scanned. Previously "marked" and "computed" "land" elements of OCEAN will retain a "marked" status at Ocean_Update re-entry, so that the user can easier keep track of previous land/ocean changes. While staying in MP, the OCEAN attributes of the elements where the cursor is located changes to "marked" (e.g. "land - marked", "site on land - marked"). + and - After MP has been entered, '+' and '-' entries function like 'l' and 'o'. s The s-key allows to return to SP mode. The mode display, however, shows that MP mode is pending, i.e. valid results will only be obtained after pressing the g-key. u and c The user should set those MRD flags TRUE where Schwiderski's tables show values that are NOT UNDERLINED. For this purpose, the cursor is moved to any place inside the corresponding 1x1 deg quadrangle of OCEAN, and the u-key is pressed. The c-key is the counterpart of the u-key. Does the model contain unconstrained 1x1 degree quadrangles with "land" AND "sea" inside ? If not, all MRD-flags should remain false. Is there very little land in an unconstrained box (small islands) ? If yes: Don't worry, the MRD effect will be insignificant. m Starts multi-point mode in rectangular areas. Initially, the rectangle contains only the cell with the curser inside. Use arrow keys to adjust its size. To accept the rectangle, press an edit function key (+ - o l), designating the desired state change. To ignore, an illegal key should be pressed (e.g. z). To shift the starting corner to the current end corner, press m again. 2.4.3.3 Convenience keys under Ocean_Update Find a color: f Permanent change of palette: p Normal cursor movement interval: n Automatic cursor increment on/off: a f The legend is displayed. The Cursor_Prompt routine scans the color tags of the labels. The user locates the cursor inside the color tag and presses the SpaceBar. If an element of OCEAN with the corresponding attribute can be found, the cursor will disappear, and the user can perform Change_Palette operations (c.f. 2.4.1.2) to find the corresponding OCEAN elements = those that are changing their color (temporary palette change). If no such element is contained in OCEAN, a beep will sound and the cursor will reappear on the OCEAN grid. p Cursor area is like under the f-key. To select a permanent palette color, c.f. Chapter 2.4.1.2 n After F1 or F2 have been used to change the cursor increment, entering the n-key will reset the cursor interval and location appropriate w.r.t. the OCEAN and screen resolution. a Normally, the keys + - o l would have to be pressed alternately with the the cursor movement keys if a long row of OCEAN elements are going to be changed/marked. For this purpose, an automatism was devised to move the curser after any of the above keys have been pressed if the cursor movement direction and increment was the same in the last two moves. Pressing the a-key toggles this option on/off. 2.4.4 UNIX X-windows version Changes affect cursor_prompt which is called under update of TOPO and OCEAN. More functionality: TOPO: Enter update mode, using the U key. You can use the mouse to locate a box. Press B to enter multi-box mode. You can use the mouse to widen/narrow the multi-box rectangle. OCEAN: You can use the mouse to locate a box. Press M to enter multi-box mode. You can use the mouse to widen/narrow the multi-box rectangle. Reduced functionality: The INKEY$() function is not available. No equivalent could yet be made. Therefore, the blinking function in Find ( x-key or f-key in OCEAN update) does not work. 2.5 Automatic processing Update the .ins-file To activate, set in namelist ¶m QA=.TRUE. or from the terminal enter site number to process followed by " A" E.g. 23 A The interactive stages of updating the TOPO array and updating the OCEAN hi-res. arrays will be skipped. Also coefficients output will be produced unconditionally. In namelist ¶m iuins=4 will redirect the site number input to the .ins-file. Then, processing can be totally automatic. However, the user is always prompted to accept/update new TOPO arrays if in namelist ¶m qptopo=.true. The totally automatic processing is recommended if large sets of sites have to be rerun for e.g. a variant of the ocean model, Green's function etc. 2.6 Saving the results Completed post-processed results are obtained after quit in Ocean_Update (1) if OLEFF has computed on the raw OCEAN grid, succeeded by SP mode editing in Ocean_Update if necessary (2) if computation was ordered from within Ocean_Update (g-key) Programmed precaution: OLMPP checks that at least one of (1) and (2) was done. If not, the loading results are certainly incomplete. The user gets a prompt, suggesting to return to the computation stage / Ocean_Update. If the user forces output of the results to file, the records will be marked "yet incomplete". OLMPP, however, cannot know whether a pending MP mode contained significant updates of OCEAN, or whether the user just played around a little before quitting the site. The output records are always appended to the file FOUT. Their format follows #BLOKQ conventions (VLBI). For the user's convenience to find the post-processed records (and delete redundant / erroneos results), an ID-label / Time Stamp option was introduced: The system prompts: Output results ((y)/N) y Enter a string to label this set (optional). Choose or combine from three forms: :text... for new label (max 48 chars.) or / " empty " or T " date & time stamp (default). Any other entry: Use previous label. MMMM DD, YYYY HH:MM _Space for text with stamp__ One could enter e.g. T elastic earth model After writing the results to file, the user has the opportunity to reprocess the currently defined OCEAN grid (e.g. redefining the Mass Redistribution flags). If the opportunity is not used, all changes to OCEAN will be lost. What you can do to close the gap ================================ Especially considering gravity effects, the coastline resolution might not be good enough even if you go to resolution stage 4. Here is a recipe what you could do in that case. (1) Consider that neither the tide nor the "normalized" Green's function vary with any considerable amount when the separation between load and field point is less than 5 km. (2) Inspect while operation on the Ocean Grid: The tide height that is interpolated at the cells closest to the site. Note down these values (ZT). Change the displayed partial tide by pressing R and then a key 1..9AB (3) Flipping the cell state to land and back to sea (keys - and + ) causes the oleff.f subroutine to print e.g. >> Using the quasi-analytical formula at near-range >> Loading : (-9.05335E-05,.0) >> Attraction [mGal]: 4.175892177173566E-02 >> Load+Newton(SeaLvl): (-2.19377E-05,.0) >> Newton at Altitude: -1.650844678730090E-04 Note down the factors for Loading and Load+Newton(SeaLvL) (GL). Remember that both the loading factor and the Newtonian factor in the near-range case represent the effect of a whole Ocean Grid cell with the gravimeter above or at its very centre. Turn the cells to land and exit the Ocean Grid. (4) That step is really hands-on. Take a map with good resolution (1:50000) Draw the ocean grid cells in the map; denote the cells that were sea and land, respectively. The missing effect is the sum of GL*ZT plus the following integral -BG*RHO*ZT*ALT* integral{ dx dy / (ALT^2 + X^2 + Y^2 )^(3/2) } separately for each partial tide 3.0 Recent updates Feb 1993: --------- Problem: System hangs after PrintScreen / Screen Capture. While system is in graphic mode, many read-keyboard-direct polls and prompts (hopefully all) have been extended to support the DOS PrintScreen service in a more safe way. Press ^P before PrintScreen, so the system receives the PrintScreen request while talking to DOS. Pressing CarriageReturn when the image is ready will take you back. Several new commands in Ocean_Update, minor importance. Check Help. Addition to Chapter 2.4.1.2, suitable for "f"-find a box: Controls have been augmented with a blink option. Boxes of the requested type (most often: where is the site ?) blink as an eye-catcher. Since blinking frequency is dependent on how fast your system reads the keyboard (INKEY$-function), the speed can be adjusted; press "+" or "-". Permanent palette change is no longer an option under the "f" command. Use "p" instead. The help text is now printed on a graphic screen. In order to output the text to a printer, press ^P, then PrintScreen key; after the image is printed, press CarriageReturn. If colors present a problem to the printer, adjust the colors of the help page by pressing "f" for foreground, "b" for background, or hit "w" to get white-on-black. October 1994: ------------- Generalize to ocean models on grids not necessarily 1x1 deg aligned on i+0.5 deg. Accept e.g. Le Provost: 0.5 x 0.5 deg aligned on i/2 deg. Speed-up interpolation with OLINTPG for the following situations: (1) The entire area is SEA. All interpolation will be on fully occupied 5x5 suspension arrays. Store a set of coefficients once for all, don't check geometry more than once. (2) For every fully occupied 5x5 suspension array re-use the set of coefficients. Geometry must be checked. There is a new ENTRY ZINTPV_RESET (l1,i1). If another interpolation package is used, a corresponding dummy subroutine must be included in order to prevent that the system goes into a trap. Use LOGICAL l1 and INTEGER i1. Use of the speed-up feature is initiated by default. A new NAMELIST variable i25U controls the feature. Adapt OLF.F to UNIX: call it olfg.f. Exception table production with deferred TOPO has the advantage that users not considering post-processing can do without TOPO. New option to save the 5x5 minutes resolution, 3x3 degrees size local TOPO grids of each site in files (*.TOP). They are preferably kept in the subdirectories pertaining to each project (GPS, AUS, ...). A record of their geographic locations and hard disk whereabouts is kept in a directory file TOPSDP/'tops.dir' The aim is to avoid multiples of files for a single site, so the files may be retrieved by site location (TOPSDP is set by DATA, variable through NAMELIST). November 1994: -------------- Interactive Update capability for TOPO array. Automatic processing (new Namelist parameters QA and IUINS). May to July 1995: ----------------- An implementation of OLMPP on an X-Windows UNIX platform has been accomplished using Caltech's pgplot package as an interface. Text display is really bad, I apologize. That's the way to install: (1) Obtain pgplot from astro.caltech.edu ftp://astro.caltech.edu/pub/pgplot and install following the instructions until "make all" but not "make clean". Suppose the library is /usr/local/pgplot/libpgplot.a (2) Obtain Olmpp.tar This archive contains updates to pgplot (hgs:pgplot) Install by copying (cp) and touching (touch) the files from hgs:pgplot/src and hgs:pgplot/drivers to /usr/local/pgplot/src and /usr/local/pgplot/drivers and "make all" again. Then "make clean". (3) Install the files in $HOME/pgplot/4mw (cp) and Fortran-compile them; collect objects in pgplot/libgrex.a (ar) Beware platform dependent date and time routines. Beware calling these routines; special compiler options ( +E1 ) may be necessary. (4) Install Oload/p and Oload/p/mpp, compile, lib and link olmpp. Observe the following files and symbolic links: ls -l /usr/local/pgplot -rw-r--r-- 1 hgs users 5207 May 19 07:32 drivers.list lrwxr-xr-x 1 hgs users 29 May 26 14:03 grex.fh -> /users/hgs/pgplot/4mw/grex.fh -rw-r--r-- 1 hgs users 66020 May 19 08:04 grfont.dat -rw-r--r-- 1 hgs users 4102 May 19 07:32 grpckg1.inc lrwxr-xr-x 1 hgs users 33 Jun 2 07:01 keys.fh -> /users/hgs/pgplot/4mw/keys.fh lrwxr-xr-x 1 hgs users 31 May 26 14:02 libgrex.a -> /users/hgs/pgplot/4mw/libgrex.a -rw-r--r-- 1 hgs users 506376 Jun 11 11:22 libpgplot.a -rwxr-xr-x 1 hgs users 385024 Jun 11 11:22 libpgplot.sl -rw-r--r-- 1 hgs users 33611 Jun 14 12:14 makefile -rw-r--r-- 1 hgs users 132042 Jun 11 11:22 pgplot.doc -rw-r--r-- 1 hgs users 129083 May 19 08:51 pgplot.html -rw-r--r-- 1 hgs users 3831 May 19 07:32 pgplot.inc -rwxr-xr-x 1 hgs users 49152 May 19 08:05 pgxwin_server -rw-r--r-- 1 hgs users 16059 May 19 07:32 rgb.txt This OLMPP manual still contains remains from the PC-version. There are only a few changes. Most things concerning files can be treated equivalently, using "/" for "\" or ignoring "D:" etc. Many character variables have been redefined, so longer file names can be specified. C.f. Chapter 2.4.4 July 1995: ---------- Fitting the Flather grid. Now one can process Flather's North Atlantic model with olfg and post-process coastal sites. All necessary information is contained in the .xtb-file. The higher resolution demanded to change the display of the MRD flags at the right margin of the OCEAN grid. Now we use = for T (true, MRD enabled) . for F (false, MRD disabled) APPENDIX A ========== olmpp.needs Shows external references fcs Compiles subroutines (instead of "make") There is a general procedure in $HOME/bin If necessary, directories may have their own fcs. fcs sub xyz where sub source file name specified without trailing ".f" xyz lib file name without leading "lib" and trailing ".a" fcm Compiles mains. Each mains directory ( examples: /m /mpp ) has its particular fcm. The call is simple: fcm main [ opt ] Nonempty option opt asks for a map file output. $HOME/bin/fcs #!/bin/sh echo Compiling $1.f with options -c -K +O2 -V ${F77OPT:=} f77 -c -K +O2 -V $F77OPT $1.f if [ x${2}x != xx -a -r $1.o ] then echo Archiving to lib$2.a ar rv lib$2.a $1.o fi pgplot/4mw/fcs compiles the MicroWay fgrex interface routines #!/bin/sh echo Compiling $1.f f77 -c +O2 +E1 +U77 -V $1.f if [ x${2}x != xx -a -r $1.o ] then echo Archiving to lib$2.a ar rv lib$2.a $1.o fi Oload/p/mpp/fcm #!/bin/sh if [ $# -eq 2 ] then fort77 -K +O2 -V +U77 +E1 +E6 ${1}.f \ -L${HOME}/Oload/p -lload \ -L${HOME}/Oload/p/gra -lgra \ -L${HOME}/Oload/p -lload \ -L${HOME}/math -lrcp \ -L${HOME}/util -lutil \ -L/usr/local/pgplot -lgrex -lpgplot \ -L/usr/lib/X11R5 -lX11 \ -o ${HOME}/bin/${1} \ -Wl,-m,-t,-v > ${1}.map echo "Load map produced" else fort77 -K +O2 -V +U77 +E1 +E6 ${1}.f \ -L${HOME}/Oload/p -lload \ -L${HOME}/Oload/p/gra -lgra \ -L${HOME}/Oload/p -lload \ -L${HOME}/math -lrcp \ -L${HOME}/util -lutil \ -L/usr/local/pgplot -lgrex -lpgplot \ -L/usr/lib/X11R5 -lX11 \ -o ${HOME}/bin/${1} fi echo "$1 with libload ready, out to ${HOME}/bin/${1}" Note that we have the graphics procedures in /gra For that reason we need the libload.a reference a second time.