Re: [ferret_users] Using xbox and zbox to find total transport

Hi Ryo,

As you suggested I tried to rectify the following script.

use "Data/Test.nc"
LIST/precision=6 ZBOXLO[GZ=v],ZBOXHI[GZ=v],Z[GZ=v],ZBOX[GZ=v]
LET NEW_EDGES=ZCAT ({0}, zboxhi[gz=v])
LIST NEW_EDGES
LET zpts = z[gz=v]
DEFINE AXIS/Z/UNITS=meters/DEPTH/BOUNDS zax = zpts, NEW_EDGES

LET newv=v[gz=zax@asn]
LIST/precision=6 ZBOXLO[GZ=newv], ZBOXHI[GZ=newv], Z[GZ=newv], ZBOX[GZ=newv]

I am also attaching the outputs with more precision, before and after correcting zaxis. Is this proper now?

DEPTH (m): -0.02 to 5500
Column 1: ZBOXLO is ZBOXLO (axis DEPTH)
Column 2: ZBOXHI is ZBOXHI (axis DEPTH)
Column 3: Z is Z (axis DEPTH)
Column 4: ZBOX is ZBOX (axis DEPTH)
               ZBOXLO   ZBOXHI      Z      ZBOX
5       / 1:    -0.02    10.04     5.01   10.060
15.1    / 2:    10.04    20.18    15.07   10.135
25.3    / 3:    20.18    30.52    25.28   10.345
35.8    / 4:    30.52    41.19    35.76   10.665
46.6    / 5:    41.19    52.30    46.61   11.110
58      / 6:    52.30    64.00    57.98   11.705
70      / 7:    64.00    76.47    70.02   12.470
82.9    / 8:    76.47    89.92    82.92   13.450
96.9    / 9:    89.92   104.62    96.92   14.700
112.3   / 10:   104.62   120.91   112.32   16.285
129.5   / 11:   120.91   139.23   129.49   18.320
149     / 12:   139.23   160.18   148.96   20.955
171.4   / 13:   160.18   184.59   171.40   24.415
197.8   / 14:   184.59   213.63   197.79   29.040
229.5   / 15:   213.63   248.97   229.48   35.335
268.5   / 16:   248.97   293.05   268.46   44.085
317.6   / 17:   293.05   349.52   317.65   56.465
381.4   / 18:   349.52   423.65   381.39   74.130
465.9   / 19:   423.65   522.61   465.91   98.960
579.3   / 20:   522.61   654.33   579.31 131.720
729.4   / 21:   654.33   823.86   729.35 169.530
918.4   / 22:   823.86 1028.76   918.37 204.900
1139.2 / 23: 1028.76 1258.86 1139.15 230.100
1378.6 / 24: 1258.86 1502.13 1378.57 243.275
1625.7 / 25: 1502.13 1750.40 1625.70 248.270
1875.1 / 26: 1750.40 2000.06 1875.11 249.655
2125    / 27: 2000.06 2250.01 2125.01 249.945
2375    / 28: 2250.01 2500.00 2375.00 249.995
2625    / 29: 2500.00 2750.00 2625.00 250.000
2875    / 30: 2750.00 3000.00 2875.00 250.000
3125    / 31: 3000.00 3250.00 3125.00 250.000
3375    / 32: 3250.00 3500.00 3375.00 250.000
3625    / 33: 3500.00 3750.00 3625.00 250.000
3875    / 34: 3750.00 4000.00 3875.00 250.000
4125    / 35: 4000.00 4250.00 4125.00 250.000
4375    / 36: 4250.00 4500.00 4375.00 250.000
4625    / 37: 4500.00 4750.00 4625.00 250.000
4875    / 38: 4750.00 5000.00 4875.00 250.000
5125    / 39: 5000.00 5250.00 5125.00 250.000
5375    / 40: 5250.00 5500.00 5375.00 250.000

------------------------------------------------------------------------------------------------------

            DEPTH (m): 0 to 5500
Column 1: ZBOXLO is ZBOXLO (axis ZAX)
Column 2: ZBOXHI is ZBOXHI (axis ZAX)
Column 3: Z is Z (axis ZAX)
Column 4: ZBOX is ZBOX (axis ZAX)
               ZBOXLO   ZBOXHI      Z      ZBOX
5       / 1:     0.00    10.04     5.01   10.040
15.1    / 2:    10.04    20.18    15.07   10.135
25.3    / 3:    20.18    30.52    25.28   10.345
35.8    / 4:    30.52    41.19    35.76   10.665
46.6    / 5:    41.19    52.30    46.61   11.110
58      / 6:    52.30    64.00    57.98   11.705
70      / 7:    64.00    76.47    70.02   12.470
82.9    / 8:    76.47    89.92    82.92   13.450
96.9    / 9:    89.92   104.62    96.92   14.700
112.3   / 10:   104.62   120.91   112.32   16.285
129.5   / 11:   120.91   139.23   129.49   18.320
149     / 12:   139.23   160.18   148.96   20.955
171.4   / 13:   160.18   184.59   171.40   24.415
197.8   / 14:   184.59   213.63   197.79   29.040
229.5   / 15:   213.63   248.97   229.48   35.335
268.5   / 16:   248.97   293.05   268.46   44.085
317.6   / 17:   293.05   349.52   317.65   56.465
381.4   / 18:   349.52   423.65   381.39   74.130
465.9   / 19:   423.65   522.61   465.91   98.960
579.3   / 20:   522.61   654.33   579.31 131.720
729.4   / 21:   654.33   823.86   729.35 169.530
918.4   / 22:   823.86 1028.76   918.37 204.900
1139.2 / 23: 1028.76 1258.86 1139.15 230.100
1378.6 / 24: 1258.86 1502.13 1378.57 243.275
1625.7 / 25: 1502.13 1750.40 1625.70 248.270
1875.1 / 26: 1750.40 2000.06 1875.11 249.655
2125    / 27: 2000.06 2250.01 2125.01 249.945
2375    / 28: 2250.01 2500.00 2375.00 249.995
2625    / 29: 2500.00 2750.00 2625.00 250.000
2875    / 30: 2750.00 3000.00 2875.00 250.000
3125    / 31: 3000.00 3250.00 3125.00 250.000
3375    / 32: 3250.00 3500.00 3375.00 250.000
3625    / 33: 3500.00 3750.00 3625.00 250.000
3875    / 34: 3750.00 4000.00 3875.00 250.000
4125    / 35: 4000.00 4250.00 4125.00 250.000
4375    / 36: 4250.00 4500.00 4375.00 250.000
4625    / 37: 4500.00 4750.00 4625.00 250.000
4875    / 38: 4750.00 5000.00 4875.00 250.000
5125    / 39: 5000.00 5250.00 5125.00 250.000
5375    / 40: 5250.00 5500.00 5375.00 250.000

I did the ncdump of my data file and I think dz information is not there in the header file, here I am also attaching the same.

dimensions:
    lon = 720 ;
    lat = 330 ;
    depth = 40 ;
variables:
    float temp(depth, lat, lon) ;
        temp:long_name = "TEMPERATURE" ;
        temp:units = "deg. C" ;
        temp:_FillValue = -9.99e+33f ;
    float salt(depth, lat, lon) ;
        salt:long_name = "SALINITY" ;
        salt:units = "g/kg" ;
        salt:_FillValue = -9.99e+33f ;
    float u(depth, lat, lon) ;
        u:long_name = "ZONAL VELOCITY" ;
        u:units = "m/sec" ;
        u:_FillValue = -9.99e+33f ;
    float v(depth, lat, lon) ;
        v:long_name = "MERIDIONAL VELOCITY" ;
        v:units = "m/sec" ;
        v:_FillValue = -9.99e+33f ;
    float w(depth, lat, lon) ;
        w:long_name = "VERTICAL VELOCITY" ;
        w:units = "m/sec" ;
        w:_FillValue = -9.99e+33f ;
    float taux(lat, lon) ;
        taux:long_name = "TAUX" ;
        taux:units = "N/m^2" ;
        taux:_FillValue = -9.99e+33f ;
    float tauy(lat, lon) ;
        tauy:long_name = "TAUY" ;
        tauy:units = "N/m^2" ;
        tauy:_FillValue = -9.99e+33f ;
    float ssh(lat, lon) ;
        ssh:long_name = "SEA SURFACE HEIGHT" ;
        ssh:units = "m" ;
        ssh:_FillValue = -9.99e+33f ;
    double lon(lon) ;
        lon:units = "degrees_east" ;
    double lat(lat) ;
        lat:units = "degrees_north" ;
    double depth(depth) ;
        depth:units = "meters" ;
        depth:positive = "down" ;

On Thu, Aug 10, 2017 at 7:47 PM, Ryo Furue <furue@xxxxxxxxxx> wrote:

Hi James,

I was trying to calculate the transport using v[x=@din,z=@din], and i compared these values with v*xbox*zbox followed by x=@sum,z=@sum]. The results I am getting is not at all comparable. I understood from the forum that depth axis with irregular dz can cause wrong result.

Will cause wrong results. I have to guess, but from your listing of ZBOX, it seems that the axis isn't properly set. For example, the first depth is 5 m. Isn't the dz of the first level 10 m? Your ZBOX says it's 10.1m.

The key is that ZBOX should be identical to the correct dz for each depth level.

Conversely, if ZBOX values are identical to the correct dz's, you are fine. (It seems you need more precision when LISTing.) An easier method to check this, is to use the ncdump command on your command line. Then, you will immediately know whether the depth axis has the correct dz information.

Then which way I should follow?.

Both ways will give you wrong results.

The "correct" way is to fix your vertical axis, if it's not properly set.

yes? define axis/z/depth/units=meters/from_data/edges zax =
{0, 10, 20, 300, . . . . , 6000}
yes? let v_fixed = v[gz=zax@ASN]

where {0, 10, 20, . . . , 6000} are the edges of the vertical grid cells. Don't use these numbers; I've just made them up. When there are N cells in the vertical, there are N+1 edges. The above DEFINE AXIS command constructs a vertical axis that consists of N points which are midpoints between the edges and adds the positions of the edges as a piece of auxiliary information to the axis. Note the EDGES qualifier.

By the way, if your goal is an integration, you don't need to compute v * xbox * zbox as an intermediate variable. The operator @DIN uses the xbox and zbox information correctly. But, a caveat is that you need to use the correct upper and lower bounds for the integration: In

yes? v_din = v[z=0:6000@DIN] ! if the total depth is 6000m.
yes? v_box = v * zbox
yes? v_box_sum = v_box[k=1:N@SUM]

v_din and v_box_sum should be the same (I think. I'm not familiar with @SUM).

Cheers,

Ryo