Re: [ferret_users] LEVITUS Ocean Data/FERRET shows differing results

[Ansley Manke <ansley.b.manke@xxxxxxxx>]

Hi Mubashar,
I thought about it some more.  In this case it is the location of missing data that's causing the difference in how these averages are computed.

Here's an example that I think shows what is happening with these averages. It comes from the fact that you're in a basin where at different XY locations there are different number of Z levels with good data.  It will be up to you to decide what calculation to use, based on what kind of average operation you want to do with this data.


If we were out in the open ocean, where there is data all the way down, you would not see the difference.

yes?  use levitus_climatology


  ! Define an XY average and a Z average variable

yes? let zave = temp[k=1:16@ave]
yes? let xyave = temp[x=330:340@ave,y=-5:5@ave]

  ! Compare the z average of the XY average, and the XY average of the Z ave.

yes? list xyave[k=1:16@ave],  zave[x=330:340@ave,y=-5:5@ave]
             DATA SET: /home/users/tmap/ferret/linux/fer_dsets/data/levitus_climatology.cdf
             LONGITUDE: 30W to 20W
             LATITUDE: 5S to 5N
             DEPTH (m): 0 to 1750
 Column  1: XYAVE[Z=@AVE] is TEMP[X=330:340@AVE,Y=-5:5@AVE]
 Column  2: ZAVE[X=@AV4,Y=@AV4] is TEMP[K=1:16@AVE]
          XYAVE   ZAVE
I / *:     7.185  7.185

Look at a slice of data in the region you're working with. At some locations the data goes deeper than others:

yes? list/y=19.5 /k=1:16/x=35e:45e temp
             VARIABLE : TEMPERATURE (DEG C)
             FILENAME : levitus_climatology.cdf
             FILEPATH : /home/users/tmap/ferret/linux/fer_dsets/data/
             SUBSET   : 10 by 16 points (LONGITUDE-DEPTH (m))
             LATITUDE : 19.5N
                 35.5E  36.5E  37.5E  38.5E  39.5E  40.5E  41.5E  42.5E  43.5E  44.5E
                 16     17     18     19     20     21     22     23     24     25
 0       /  1:   ....   ....  28.90  29.12  29.25  29.28   ....   ....   ....   ....
 10      /  2:   ....   ....  28.76  28.98  29.12  29.15   ....   ....   ....   ....
 20      /  3:   ....   ....  28.59  28.78  28.85  28.78   ....   ....   ....   ....
 30      /  4:   ....   ....  28.38  28.54  28.55  28.38   ....   ....   ....   ....
 50      /  5:   ....   ....   ....  27.80  27.70  27.39   ....   ....   ....   ....
 75      /  6:   ....   ....   ....  26.10  26.15  26.08   ....   ....   ....   ....
 100     /  7:   ....   ....   ....  24.50  24.49   ....   ....   ....   ....   ....
 150     /  8:   ....   ....   ....  22.68  22.60   ....   ....   ....   ....   ....
 200     /  9:   ....   ....   ....  21.98  21.98   ....   ....   ....   ....   ....
 300     / 10:   ....   ....   ....  21.67  21.54   ....   ....   ....   ....   ....
 400     / 11:   ....   ....   ....  21.57  21.42   ....   ....   ....   ....   ....
 600     / 12:   ....   ....   ....  20.71  20.58   ....   ....   ....   ....   ....
 800     / 13:   ....   ....   ....  21.65   ....   ....   ....   ....   ....   ....
 1000    / 14:   ....   ....   ....  21.67   ....   ....   ....   ....   ....   ....
 1200    / 15:   ....   ....   ....   ....   ....   ....   ....   ....   ....   ....
 1500    / 16:   ....   ....   ....   ....   ....   ....   ....   ....   ....   ....


I'll define a simpler variable with more extreme variation so we can really see what's going on.

  ! define a synthetic variable in X and Z

yes? let var = {41,42,43,44,45,,32,33,34,35,,22,23,24,,,,13,14,}
yes? def axis/z/depth zax =  {5,15,25,50}
yes? def axis/x=31:39:2/units=degrees_east  xax
yes? let xzvar = x[gx=xax]+z[gz=zax]
yes? let var_xz = reshape(var,xzvar)

yes? list var_xz
             VARIABLE : RESHAPE(VAR,XZVAR)
             SUBSET   : 5 by 4 points (LONGITUDE-Z)
             31E    33E    35E    37E    39E  
              1      2      3      4      5
 5    / 1:  41.00  42.00  43.00  44.00  45.00
 15   / 2:   ....  32.00  33.00  34.00  35.00
 25   / 3:   ....  22.00  23.00  24.00   ....
 50   / 4:   ....   ....  13.00  14.00   ....

  ! The average over x and the average over Z

yes? let xave =  var_xz[x=30:40@ave]
yes? let zave =  var_xz[k=1:4@ave]

  ! Here is the depth average, each column averaged. The locations with more
  ! depth data have lower values

yes? list/order=x zave
             VARIABLE : VAR_XZ[K=1:4@AVE]
             SUBSET   : 5 points (LONGITUDE)
             Z        : 0 to 62.5
           31E    33E    35E    37E    39E  
            1      2      3      4      5
          41.00  30.00  23.80  24.80  40.00

  ! Here is the average across each row.

yes? list xave
             VARIABLE : VAR_XZ[X=30:40@AVE]
             SUBSET   : 4 points (Z)
             LONGITUDE: 30E to 40E
 5    / 1:  43.00
 15   / 2:  33.50
 25   / 3:  23.00
 50   / 4:  13.50

  ! The average of the column-average is not the same as the average of the row-average.

yes? list/noheader zave[x=30:40@ave], xave[k=1:4@ave]
I / *:     31.92  24.08

On 3/25/2013 11:41 AM, Ansley Manke wrote:

Hi,
This is discussed in the documentation. All Ferret calculations take into account data marked as missing.

I think the differences you are seeing has to do with the order of operations, and the weighting that the average transformation does. Averages use cell size to compute a weighted average, and also apply a correction for the variation in the size of longitude/latitude grid boxes depending on the location in latitude on the Earth.  If it were an unweighted average then the order of operations would not matter.

You could do an integral as Olivier suggests, but I don't think that's necessary.  Do it like this, and Ferret will weight the calculation by the size of the grid boxes or parts of grid boxes in the region you have specified, where there is valid data:

     let redSea = t_an[k=1:16@ave,x=30E:45E@ave,y=10N:30N@ave]

Here are some of the relevant pieces of documentation:

http://ferret.pmel.noaa.gov/Ferret/documentation/users-guide/variables-xpressions/XPRESSIONS#_VPINDEXENTRY_439

@AVE—average

The transformation @AVE computes the average weighted by grid box size—a single number representing the average of the variable between two endpoints.

If @AVE is specified simultaneously on multiple axes the calculation will be performed as a multiple integration rather than as sequential single integrations. The output will document this fact by showing @AV4 or "XY ave" as the transformation.

See the General Information on transformations for important details about this transformation. In particular, note the discussion about specifying the averaging interval in world coordinates, e.g. longitude, meters, time as in var[x=3.4:4.6@AVE] versus specifying the interval using indices, as in var[I=4:12@AVE]. When the interval is expressed in world coordinates, the weighting is done using partial grid boxes at the edges of the interval. If the interval is expressed using indices, the entire grid cells contribute to the weights.

That link to the General Information about transformations also contains important information.

How does Ferret compute averages and integrals on the surface of the surface of the earth?

Ansley


On 3/25/2013 9:53 AM, Star Physicist wrote:

Do you think if my NETCDF file has missing value flag, FERRET considers it in calculation? Could you please slightly elaborat ewhy do we need RedSea_Mask and how it is working to compute volume? I am pasting some part of ncdump -h of my data as follows:

ncdump -h tanom_1955-2004_missing_value.ncA/TANOM/ftp.nodc.noaa.gov/pub/WOD/SELECT/CDOME
netcdf tanom_1955-2004_missing_value {
dimensions:
    gsize = 2 ;
    lon = 360 ;
    nb2 = 2 ;
    lat = 180 ;
    depth = 16 ;
    time = UNLIMITED ; // (200 currently)
variables:
    float lon(lon) ;
        lon:standard_name = "longitude" ;
        lon:long_name = "longitude" ;
        lon:units = "degrees_east" ;
        lon:axis = "X" ;
        lon:bounds = "lon_bnds" ;
    float lon_bnds(lon, nb2) ;
    float lat(lat) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "latitude" ;
        lat:units = "degrees_north" ;
        lat:axis = "Y" ;
        lat:bounds = "lat_bnds" ;
    float lat_bnds(lat, nb2) ;
    float depth(depth) ;
        depth:units = "meters" ;
        depth:axis = "Z" ;
        depth:bounds = "depth_bnds" ;
    float depth_bnds(depth, nb2) ;
    double time(time) ;
        time:standard_name = "time" ;
        time:units = "months since 1955-01-01 00:00:00" ;
        time:calendar = "standard" ;
    float climatology_bounds(time, gsize) ;
        climatology_bounds:comment = "This variable defines the bounds of the climatological time period for each time" ;
    float t_an(time, depth, lat, lon) ;
        t_an:long_name = "Objectively Analyzed Mean" ;
        t_an:units = "degrees_celsius" ;
        t_an:comment = "Objectively analyzed climatologies are the objectively interpolated mean fields for an oceanographic variable at standard depth levels for the World Ocean." ;
        t_an:cell_methods = "area:mean depth:mean time:mean" ;
        t_an:missing_value = 9.96921e+36f ;

Thanks,

Mubashar

---

From: olivier.marti@xxxxxxxxxxxx
Subject: Re: [ferret_users] LEVITUS Ocean Data/FERRET shows differing results
Date: Sun, 24 Mar 2013 18:44:03 +0100
To: ferret_users@xxxxxxxx; couragetoknowgc@xxxxxxxxxxx

Le 24 mars 2013 à 18:24, Star Physicist a écrit :

Which way I should proceed then? How could I overcome this issue? Actually there is a function in CDO for vertical averaging over all grid points called "cdo vertmean ifile ofile" but I don't know whether it can handle bottom shallower sea level. Any suggestion?

You should compute a full 3D integral, not 1D/2D ou 2D/1D :

USE tanom_1955-2004_missing_value.nc

! Compute a Red Sea mask. Data dependant, what do you have for missing value ?

LET RedSea_mask = IF t_an GT -100 THEN 1 ! Use an ad hoc criterion, depending of your data

! Integrate temperature

LET RedSea_Temp = t_an[K=1:16@DIN,x=30E:45E@DIN,y=10N:30N@DIN]

! Integrate volume

LET RedSea_Vol = RedSea_mask[K=1:16@DIN,x=30E:45E@DIN,y=10N:30N@DIN]

LET anom = RedSea_Temp / RedSea_Vol 

-- 

Dr. Olivier Marti  - mailto:olivier.marti@xxxxxxxxxxxx

Laboratoire des Sciences du Climat et de l'Environnement

Institut Pierre Simon Laplace

Laboratoire CEA-CNRS-UVSQ  -  UMR 8212

Tel : +33 1 69 08 77 27 - Mob : +33 6 45 36 43 74

Fax : +33 1 69 08 30 73

<http://www.lsce.ipsl.fr/Pisp/56/olivier.marti.html>