Home >> Forums >> Going from z_w coordinate to depth in meters

Going from z_w coordinate to depth in meters


In Ichthyop there's a function in Dataset that goes from z_rho coordinates to depth in meters:


I need to get the depth associated to a certain z_w coordinate. As there's a double[][][] called
z_w_tp0 that says it contains the depth at w point, taking account of free surface elevation, I wrote the following function (I'm talking about ROMS3D):


        public double zw2depth(double x, double y, double z) {

            final double kz = Math.max(0.d, Math.min(z, (double) nz ));        
            final int i = (int) Math.floor(x);
            final int j = (int) Math.floor(y);
            final int k = (int) Math.floor(kz);
            double depth = 0.d;
            int layerAvgPoints=2;
            final double dx = x - (double) i;
            final double dy = y - (double) j;
            final double dz = kz - (double) k;
            double co = 0.d;
            double z_w;
            for (int ii = 0; ii < 2; ii++) {
                for (int jj = 0; jj < 2; jj++) {
                    for (int kk = 0; kk < layerAvgPoints; kk++) {
                        co = Math.abs((1.d - (double) ii - dx)
                                * (1.d - (double) jj - dy)
                                * (1.d - (double) kk - dz));
                        if (isInWater(i + ii, j + jj)) {
                            z_w = z_w_tp0[k + kk][j + jj][i + ii];
                            depth += co * z_w;
            return depth;

At this point, my questions are two:

- Does this code makes sense? Anyone sees any obvious mistake or something?
- Why all this spatial averaging? I observed that this is a technique widely used at several places in the code, but I don't really know the reason for that.

Thank you in advance,
Ignacio Vidal

Ichthyop version: 
Ichthyop 3.1
Hydrodynamic dataset: 
Java version: 
Java SE 8
Operating system: 

Hi Ignacio,

It looks ok to me, you adapted the z2depth function to the z_w, and you added an extra condition at the surface I guess because you got an error message at some point ?

As for why do we interpolate linearly ? Just te refine the results at sub-grid scale and avoid discontinuities at cell edges. We could for sure adopt a different approach : same z_w value anywhere in a given (i, j) cell. But in order to avoid discontinuities between cells it is smoother to consider that a value at given (x, y) point is an interpolation of the 4 surrounding cells (or 8 if you take into account the vertical dimension). And we went for trilinear interpolation because it is classical but we could use other interpolation method, for instance Inverse Distance Weight. I tried that in the new HYCOM plugin that will be released soon with a parameter that let you choose how many surrounding cell points you want to use in the interpolation (8, 64, etc.)

Cheers, Philippe