|
Dear coleagues, I'm trying to write a ferret script to calculate the variation of the volume-averaged kinetic energy. My objective is to analyze some pom2k runs. The variables saved in the netcdf files are:+++++++++++++++++++++++++++++++++++++++++++++++++++++++ currently SET data sets: 1> ./r1.nc (default) name title I J K L ZZ sigma of cell centre ... ... 1:10 ... DX grid increment in x 1:175 1:150 ... ... DY grid increment in y 1:175 1:150 ... ... EAST_U easting of u-points 1:175 1:150 ... ... EAST_V easting of v-points 1:175 1:150 ... ... EAST_E easting of elevation points 1:175 1:150 ... ... EAST_C easting of cell corners 1:175 1:150 ... ... NORTH_U northing of u-points 1:175 1:150 ... ... NORTH_V northing of v-points 1:175 1:150 ... ... NORTH_E northing of elevation points 1:175 1:150 ... ... NORTH_C northing of cell corners 1:175 1:150 ... ... ROT Rotation angle of x-axis wrt. e 1:175 1:150 ... ... H undisturbed water depth 1:175 1:150 ... ... FSM free surface mask 1:175 1:150 ... ... DUM u-velocity mask 1:175 1:150 ... ... DVM v-velocity mask 1:175 1:150 ... ... RMEAN horizontally-averaged rho 1:175 1:150 1:10 ... UAB depth-averaged u 1:175 1:150 ... 1:17 VAB depth-averaged v 1:175 1:150 ... 1:17 ELB surface elevation 1:175 1:150 ... 1:17 U x-velocity 1:175 1:150 1:10 1:17 V y-velocity 1:175 1:150 1:10 1:17 W z-velocity 1:175 1:150 1:10 1:17 T potential temperature 1:175 1:150 1:10 1:17 S salinity x rho / rhoref 1:175 1:150 1:10 1:17 RHO (density-1000)/rhoref 1:175 1:150 1:10 1:17 +++++++++++++++++++++++++++++++++++++++++++++++++++++++ Follows my script: use pom2k let area=dx*dy*fsm let htot=h+elb let volcol=area*htot let vel=(u^2+v^2)^.5 let dens=(rho*1000+1000) let velrho=vel^2*dens/2 <=========== at this point I should multiply by cell volume (each one along the vertical)
let ek=velrho[x=@sum,y=@sum,z=@sum] let mke=ek/volcol plot mke
The issue is: I have no idea how to calculate the volume of each cell, because I don't know the height of each cell. The only vertical variable saved by POM is zz, the mid point between two consecutive sigma levels. I was wondering if it would be possible to set a fixed height through an axis replacement, so I did: define axis/z=0:5800:50/unit=meters/depth zdepth let v_z=zaxreplace(vel,h*zz*-1,z[gz=zdepth]) let rho_z=zaxreplace(rho,h*zz*-1,z[gz=zdepth]) Now the height is 50 m, but how could I compute the operation for each cell along the vertical axis in a varying bathymetry??? I would appreciate your
comments. If you think this script is cumbersome, suggestions would be welcomed! Thanks in advance.
Best regards, Alex
--
|