Objective Analysis Package

The objective analysis (OA) package can be used to map observations to the application grid using space-time correlation functions. The source code can be obtained as tar file from our web site.


This OA package is derived from an earlier code that I wrote in 1993, when I was at Harvard University. The basic reference for the OA equations used here is described by Carter and Robinson (1987). Comprehensive description of this methodology can also be found in Gandin (1963), Bretherton et al. (1976), Lorenc (1978), McWilliams et al. (1986), Daley (1991), among others.


This package uses a local solution for the OA equations. That is, only NNCE influential observations are considered at each mapped grid point. This method is practical because it avoids inverting large matrices when the number of observations NOBS is too large. Observations that are too far apart in space and time from the mapped point contribute very little to its estimate (correlations are very small); using all observations is too much work for such a little gain. We also have a global version of the package which uses all the observations. Please, contact me if you are interested in the global OA package.


The OA package can be used to prepare initial, climatology, and forcing conditions for our suite of hierarchy models. It also can be used to prepare fields for data assimilation. Currently, it processes the following fields:



The OA package also produces the mapping estimation error covariance and mean values for each of the above fields which are used during data assimilation.

There are several options for the computation of dynamic height. This is only possible when pressure (or depth), temperature, and salinity are available at input. These options refer to how the vertical integration of the specific volume anomaly is carried out and how the level of no-motion (ZREF) is selected:


  1. Integrate specific volume anomaly from ZREF to the surface.
  2. Integrate specific volume anomaly from ZREF to the bottom.
  3. Integrate in both directions, from ZREF to the surface and from ZREF to the bottom.
  4. Integrate specific volume anomaly from the local bottom of the hydrographic profile to the surface.
  5. Integrate specific volume anomaly from the surface to the local bottom of the hydrographic profile.
  6. Same as 3, but if the local depth is shallower than ZREF, only integrate from the local depth to the surface.

The OA mapping is carried out over the field anomaly (departure from the mean). There are several strategies on how to compute the mean field. The mean field can be computed from observations and/or climatology data. The mean can be local (a mean value at each grid point) or global (a single value for the entire grid). If the option is local, an extra OA is carried to compute the mean field but with larger space and time correlations. See ICLIMA parameter in input script oa.in.


For very sparse observation data, it is recommended to use both observation and climatology data when computing mean fields. For example, in the ocean, climatology data is a more oceanographically correct field than any extrapolation to places of no data. In dynamically active regions, the melding of observations and climatology data to compute mean fields is like assuming non-stationary statistics for the mean.


A detailed description of each of the OA parameters is provided at the bottom of input script file oa.in.

Input Data

The input data to the OA package is in our own hydrographic and atmospheric (HDAT) data format. Currently, the HDAT type file is the only data format supported in the OA package. Several codes are provided to convert to HDAT format.


I created the HDAT format in 1992, when I was at Harvard, for the same purpose. Before creating a new format, I examined several data formats available but none of them was simple, self-described, and flexible enough to be used in the OA package. The NetCDF type file was considered but discarded since in situ ocean and atmospheric data has usually a different number of observation points at each station and NetCDF does not allowed more than one unlimited dimension. Any atempt to interpolate data to regular vertical or horizontal grids to facilitate data format is against of the purpose of the OA. The OA is used to interpolate scattered raw data to application grids!


The HDAT file is an ASCII file and consists of two parts:

  1. a global header block, and
  2. a data block.

The global header contains several lines of text describing the data. The number of lines of text in the global header is unlimited. The only requirement is that the last line of text has the word END. The global header may have following information:


The data block contains the station header and data. It has the following variables:
    NHVAR           number of variables in station (integer)
    NHPTS           number of observation points (integer)
    CASTID          station identification (integer)
    HLNG            station longitude (real)
    HLAT            station latitude (real)
    >HDPTH           station elevation/depth (real)
    HTIME           time of observation (real)
    HSCLE(1:NHVAR)  unpacking data scales (real array)

    HFLAG           station flag (integer)
    HTYPE           instrument type and variables order (string)

    IHDAT(1:NHPTS)  packed data for each variable (integer array)
   

The file readhydro.F is used to read and decode a HDAT file. This routine has two entries:


The entry read_header reads and decodes the global header whereas the entry read_hydro reads and unpacks one station data per call. The station information and data is stored in common blocks in the include file hydro.h.
The HTYPE character string is decoded to the following variables:
     HTOOL          instrument (uppercase string). Possible values:

                     CTD     conductivity-temperature-depth station
                     FRC     forcing station
                     XBT     expendable bathythermograph station
                     XBTS    salted XBT station
                     XCTD    expendable CTD station

     HFLDS(1:NHVAR) variable descriptor (string array).  Possible
                     values:

                     SST     sea surface temperature
                     cf      cloud fraction
                     dQdT    net heat flux sensitivity to SST
                     er      evaporation rate
                     jd      modified Julian day
                     hF      surface net heat flux
                     lR      incoming longwave radiation
                     pair    surface air pressure
                     pr      precipitation rate
                     rh      relative humidity
                     s       salinity
                     sR      solar shortwave radiation
                     t       in situ temperature
                     tair    surface air temperature
                     tdew    surface dew-point temperature
                     us      surface U-momentum stress component
                     uw      surface U-wind component
                     vs      surface V-momentum stress component
                     vw      surface V-wind component
                     yd      year day
                     wF      surface freshwater flux, E-P
                     z       elevation/depth
   
The station header is read using the following Fortran statements:
      read (iunit,*) nhvar, nhpts, castid, hlng, hlat, hdpth, htime,
     &              (hscle(m), m = 1, nhvar)
      read (iunit,'(1x,i2,1x,a)') hflag, htype
   
The station data is read and unpack using the following Fortran statements:
      do m = 1, nhvar
        read (iunit,*) (ihdat(n), n = 1, nhpts)
        do n = 1, nhpts
          hdat(n,m) = FLOAT(ihdat(n)) * hscle(m)
        enddo
      enddo
   
where hdat is a real matrix containing the unpacked data for the current station. For more details see file readhydro.F. All the variables are documented in file hydro.h.

The following example shows an HDAT file for Levitus Annual mean climatology. Here, the number of points corresponds to number of depth points. For brevity, only three stations are shown:

title = Levitus One-Degree Annual Climatology, 1994
stations = 42164
lng_min = -179.5000
lng_max = 179.5000
del_lng = 1.0000
lat_min = -77.5000
lat_max = 89.5000
del_lat = 1.0000
format = ascii, record interleaving
type = CTD
fields = depth, temperature, salinity
units = meter, Celsius, PSU
creation_date = Wednesday - February 12, 1997 - 2:02:17 PM
END
 3     6       1  173.5000 -77.5000    75.0     0.0000 1.00E-01 1.00E-03 1.00E-03
  0 'CTD: z t s'
     0   100   200   300   500   750
  -451  -458  -556  -684  -964 -1191
 34420 34431 34456 34461 34486 34534
 3    10       2  174.5000 -77.5000   200.0     0.0000 1.00E-01 1.00E-03 1.00E-03
  0 'CTD: z t s'
     0   100   200   300   500   750  1000  1250  1500  2000
  -463  -468  -560  -682  -959 -1185 -1247 -1308 -1348 -1514
 34417 34425 34444 34448 34472 34520 34556 34575 34604 34630
 3    12       3  175.5000 -77.5000   300.0     0.0000 1.00E-01 1.00E-03 1.00E-03
  0 'CTD: z t s'
     0   100   200   300   500   750  1000  1250  1500  2000
  2500  3000
  -478  -482  -567  -684  -957 -1182 -1239 -1293 -1325 -1481
 -1659 -1818
 34412 34418 34432 34434 34459 34507 34543 34563 34591 34617
 34643 34668

The following example shows an HDAT file for UWM/COADS monthly climatology. Here, the number of points corresponds to number of time (monthly) records. For brevity, only two stations are shown:

title = UWM/COADS Global Monthly Climatology, 1945-89
stations = 42164
lng_min = -179.5
lng_max = 179.5
del_lon = 1.0000
lat_min = -77.5
lat_max = 89.5
del_lat = 1.0000
str_time = 15.000
end_time = 345.000
format = ascii, record interleaving
type = FRC, forcing fields
fields_01 = time (day, 360 days year cycle)
fields_02 = constrained incoming net heat flux (Watts meter-2)
fields_03 = incoming shortwave radiation (Watts meter-2)
fields_04 = constrained freshwater flux, E-P (centimeter day-1)
fields_05 = sea surface temperature, SST (Celsius)
fields_06 = net heat flux sensitivity to SST, dQdSST (Watts meter-2 K-1)
fields_07 = corrected surface zonal stress (Newton meter-2)
fields_08 = corrected surface meridional stress (Newton meter-2)
creation_date = Saturday - April 5, 1997 - 2:9:57.7413 PM
END
 8    12       1  173.5000 -77.5000     0.0     0.0000 1.00e+00 1.00e-01 1.00e-01 1.00e-03 1.00e-03 1.00e-01 1.00e-03 1.00e-03
  0 'FRC: yd hF sR wF SST dQdT us vs'
    15    45    75   105   135   165   195   225   255   285
   315   345
   766   -72 -1240  -543  -874 -1202  -911  -907  -761  -173
    57   967
  1664   937   306    33     0     0     0     4   169   819
  1516  1925
  -356  -342  -485  -299   -22   136   -85 -1734  -131 -1154
   -91  -374
   -80  -360 -1080 -1630 -2150 -1490  -710  -710  -710  -710
 -1410  -560
  -353  -346  -395  -268  -380  -385  -382  -382  -382  -382
  -352  -349
    -5   -13     7   -42     2    57    -8    -8    -8    -8
     2     5
    33    57    70    27    87    81    51    51    51    51
    65    48
 8    12       2  174.5000 -77.5000     0.0     0.0000 1.00e+00 1.00e-01 1.00e-01 1.00e-03 1.00e-03 1.00e-01 1.00e-03 1.00e-03
  0 'FRC: yd hF sR wF SST dQdT us vs'
    15    45    75   105   135   165   195   225   255   285
   315   345
   765   -67 -1237  -468  -870 -1199  -907  -903  -757  -170
    54   965
  1658   935   305    33     0     0     0     4   169   818
  1513  1918
  -369  -345  -490  -314   -22   121   -94 -1697  -139 -1166
   -95  -383
   -90  -360 -1080 -1570 -2150 -1490  -710  -710  -710  -710
 -1410  -560
  -353  -346  -395  -357  -381  -385  -382  -382  -382  -382
  -351  -350
    -5   -13     6    18     2    57    -8    -8    -8    -8
     2     4
    33    57    70    -7    86    81    51    51    51    51
    66    48

The input grid positions for the mapping are read from GRID NetCDF file. See grid generation tools for more information.

Output Data

The output mapped fields and their associated mean and errors fields are written out into a NetCDF file. The name of the NetCDF variables are:

Table 1
Observation Estimate Mean Variance
In situ temperature temp temp_mean temp_err
Potential temperature temp temp_mean temp_err
Surface air temperature Tair Tair_mean Tair_err
Surface dew-point temperature Tdew Tdew_mean Tdew_err
Salinity salt salt_mean salt_err
In situ density anomaly dena dena_mean dena_err
Sigma-t sigmat sigmat_mean sigmat_err
Sound speed sndsp sndsp_mean sndsp_err
Dynamic height dynht dynht_mean dynht_err
Surface air pressure Pair Pair_mean Pair_err
Sea surface height SSH SSH_mean SSH_err
Surface net heat flux shflux shflux_mean shflux_err
Surface freshwater flux (E-P) swflux swflux_mean swflux_err
Precipitation rate precip precip_mean precip_err
Evaporation rate evapor evapor_mean evapor_err
Solar shortwave radiation swrad swrad_mean swrad_err
Longwave radiation lwrad lwrad_mean lwrad_err
Relative humidity humidity humidity_mean humidity_err
Cloud fraction cloud cloud_mean cloud_err
Surface U-wind component Uwind Uwind_mean Uwind_err
Surface V-wind component Vwind Vwind_mean vwind_err
Surface U-wind stress component sustr sustr_mean sustr_err
Surface V-wind stress component svstr svstr_mean svstr_err
Sea surface temperature SST SST_mean SST_err
Net heat flux sensitivity to SST dQdSST dQdSST_mean dQdSST_err
Sea surface salinity SSS SSS_mean SSS_err


All current OA fields are mapped at the RHO-positions (center of grid box) of a staggered C-grid. The surface wind and wind stress are also at RHO-points to facilitate averaging to U- and V-points when rotating to the appropriate curvilinear grid angle.

Post-processing

ROMS initial and climatology conditions NetCDF files are generated from output OA NetCDF file by the initialization and climatology package. In the same fashion, the forcing package can be used to generate ROMS forcing NetCDF from the output OA file.

References
  1. Bretherton et al., 1976: A technique for objective analysis and design of oceanographic experiments applied to MODE-73, Deep Sea Res., 23, 559-582.


  2. Carter and A.R. Robinson, 1987: analysis models for the estimation of oceanic fields, J. Atms. and Ocean. Techn., 4, 49-74.


  3. Daley, R., 1991: Atmospheric data anlysis, Cambridge University Press, Cambridge, 457 pp.


  4. Gandin, L.S., 1963: The objective analysis of meteorological fields. Leningrad, Hydrometeorological Publ. House. English translation, Israel Program for Scientific Translations, Jerusalem, 1965, 242 pp.


  5. Lorenc, A.C., 1986: Analysis methods for numerical weather prediction, Quart. J. R. Met. Soc., 112, 1177-1194.


  6. McWilliams, et al., 1986: An objective analysis of the POLYMODE Local Dynamics Experiment. Part I: general formalism and statistical model selection, J. Phys. Oceanogr., 16, 483-504.






IMCS, Ocean Modeling Group, by Webmaster.