[Thread Prev][Thread Next][Index]

Re: [ferret_users] Problem to plot grid


I've tried that but still could not get anywhere.
I've actually added one point at each extreme, changing my grid from 150x300 to 152x302.
However, I am still getting the same problem with the last cell just being half plotted, as you can see below. The plot on the left uses the original lon_rho and lat_rho and the one on the right uses the edited lon_rho2 and lat_rho2. (I've created an over complicated script to generate these figures... for that reason I will just add it to the end of the email). when I tried to plot it just until the second to last cell (in this case i=151), then that cell was also half plotted.

Best regards,

Inline images 2

use grids_example.nc


def ax/x=0:151:1 xax1
def ax/y=0:301:1 yax1

!Generating new first and last cells
let lon_rhoa = lon_rho[gx=xax1,gy=yax1]
let lon_rhob = if i lt 152 then lon_rhoa else (lon_rhoa[i=150]+(lon_rhoa[i=150]-lon_rhoa[i=149]))
let lon_rhoc = if j lt 302 then lon_rhob else (lon_rhob[j=300]+(lon_rhob[j=300]-lon_rhob[j=299]))
let lon_rhod = if i gt 1 then lon_rhoc else (lon_rhoc[i=2]+(lon_rhoc[i=2]-lon_rhoc[i=3]))
let lon_rhoe = if j gt 1 then lon_rhod else (lon_rhod[j=2]+(lon_rhod[j=2]-lon_rhod[j=3]))

let lat_rhoa = lat_rho[gx=xax1,gy=yax1]
let lat_rhob = if i lt 152 then lat_rhoa else (lat_rhoa[i=150]+(lat_rhoa[i=150]-lat_rhoa[i=149]))
let lat_rhoc = if j lt 302 then lat_rhob else (lat_rhob[j=300]+(lat_rhob[j=300]-lat_rhob[j=299]))
let lat_rhod = if i gt 1 then lat_rhoc else (lat_rhoc[i=2]+(lat_rhoc[i=2]-lat_rhoc[i=3]))
let lat_rhoe = if j gt 1 then lat_rhod else (lat_rhod[j=2]+(lat_rhod[j=2]-lat_rhod[j=3]))

!Calculating limits between grid cells
let lon_rho2 = ((lon_rhoe[i=@shf:0,j=@shf:0]+lon_rhoe[i=@shf:1,j=@shf:0]+lon_rhoe[i=@shf:0,j=@shf:1]+lon_rhoe[i=@shf:1,j=@shf:1])/4)
let lat_rho2 = ((lat_rhoe[i=@shf:0,j=@shf:0]+lat_rhoe[i=@shf:1,j=@shf:0]+lat_rhoe[i=@shf:0,j=@shf:1]+lat_rhoe[i=@shf:1,j=@shf:1])/4)


set win/asp=0.5 2
let diff_rho = lon_rho+lat_rho
let diff_rhoe = lon_rhoe+lat_rhoe

def symb hlim -74.8:-74.55
def symb vlim 36.58:36.7
def symb lev (-38.2,-37.9,0.01)

can mode verify
set view left
shade/i=1:150/j=1:300/lev=($lev)/hlim=($hlim)/vlim=($vlim)/nolab diff_rho,lon_rho,lat_rho
rep/i=1:150 (plot/ov/vs/dash/line/col=1/nolab lon_rho,lat_rho)
rep/j=1:300 (plot/ov/vs/dash/line/col=1/nolab lon_rho,lat_rho)
rep/i=148:150 (rep/j=1:3 (plot/size=0.2/col=1/ov/vs/nolab/symb=2 {`lon_rho`},{`lat_rho`}))
label/nouser `($ppl$xlen)/2` `($ppl$ylen)*1.1` 0 0 0.2 Original variables

set view right
shade/i=1:152/j=1:302/lev=($lev)/hlim=($hlim)/vlim=($vlim)/nolab diff_rhoe,lon_rhoe,lat_rhoe
rep/i=1:151 (plot/ov/vs/dash/line/col=1/nolab lon_rho2,lat_rho2)
rep/j=1:300 (plot/ov/vs/dash/line/col=1/nolab lon_rho2,lat_rho2)
rep/i=148:150 (rep/j=1:3 (plot/size=0.2/col=1/ov/vs/nolab/symb=2 {`lon_rho2`},{`lat_rho2`}))
label/nouser `($ppl$xlen)/2` `($ppl$ylen)*1.1` 0 0 0.2 Edited variables

set mode verify

On 10 September 2015 at 21:23, Russ Fiedler <russell.fiedler@xxxxxxxx> wrote:

Hi Andre,

Yes, you will need to construct the corners of each cell somehow. This will allow Ferret to extrapolate correctly out to the edges.
If you have the corners saved, good. Then all you need to do is

Assuming your grid is 150x300 the corner arrays will be 151x301

shade var,lon_corner,lat_corner
repeat/i=1:151 plot/ov/vs/dash lon_corner,lat_corner
repeat/j=1:301 plot/ov/vs/dash lon_corner,lat_corner

If you don't have the corners saved I would write a quick program (not Ferret) to calculate them and save them to a file which can be read in.
You can do this in Ferret but it's quite fiddly. You will also need to interpolate in 2D and not 1D as you have done below.



On 11/09/15 06:22, Andre Paim wrote:
Dear Russ and Ansley, thanks a lot for your answers.

I  think I wasn't as thorough as I should have when explaining my issue. Attached is a netcdf with the the lon and lat variables.
This is the code I used to create the first figure in the previous email:


let lon_hyc2 = (lon_hyc+lon_hyc[i=@shf:1])/2
let lat_hyc2 = (lat_hyc+lat_hyc[j=@shf:1])/2


set win/asp=1 1
shade/title="Hycom Grid"/nokey/x=-75.8:-73.8/y=38.0:40.0 surf_hyc
go land_detail

rep/i=1:126 (plot/ov/vs/line/dash/col=1/nolab lon_hyc2,lat_hyc2)
rep/j=1:87  (plot/ov/vs/line/dash/col=1/nolab lon_hyc2,lat_hyc2)


The /dash in here is just because this is not my main grid. What I'm trying to accomplish here is to just show the differences between these two grids (one from the HYCOM available at hycom.org, and the other one is my grid) at my boundaries and in my area of interest.

So, the HYCOM grid works fine. However when I try to do the same to my grid (using the code below), I get the second figure


let lon_rho2 = (lon_rho+lon_rho[i=@shf:-1])/2
let lat_rho2 = (lat_rho+lat_rho[j=@shf:1])/2


set win/asp=1 2
shade/lev=(-38.5,-37.4,0.02)/hlim=-75.2:-74.4/vlim=36.4:37.1/nolab lon_rho+lat_rho,lon_rho,lat_rho

rep/i=2:150 (plot/ov/vs/line/col=7/nolab lon_rho2,lat_rho2)
rep/j=1:299 (plot/ov/vs/line/col=7/nolab lon_rho2,lat_rho2)


Inline images 2

This second figure is actually a little different from the original; I figured a mistake a had done in my code. However, the offset is still there. I've plotted using /dash just to make it easier to notice the offset.
Also, it seems that I'm not showing the first row and last column. And they seem to have just half of the size of the other cells.

Russ, I think the link you suggested me kind of cover this last issue, but I didn't understand how to use that. Should I write:

shade/lev=(-38.5,-37.4,0.02)/hlim=-75.2:-74.4/vlim=36.4:37.1/nolab lon_rho+lat_rho,lon_rho[i=1:151,j=0:300],lat_rho[i=1:151,j=0:300]

I doesn't seem to make any difference.

Any suggestions?

Thanks a lot!

On 9 September 2015 at 22:44, Russ Fiedler <russell.fiedler@xxxxxxxx> wrote:

Hi Ansley,

If you look at your plot you see lines at approximately half the cell heights. This is an artifact due to joining the wrapped endpoints (380E ->20E). Missing values need to be introduced at the ends. The vertical lines can be plotted by plotting the transposes of lon2 and lat2

yes? use coads_climatology

! Cell centres
yes? let/title="longitude"/units=degrees_east lon = x[gx=sst] + 0*y[gy=sst]
yes? let/title="latitude"/units=degrees_north  lat = 0*x[gx=sst] + y[gy=sst]
yes? shade/l=1/hlim=180:240/vlim=0:40 sst, lon, lat

! Define cell NE points
yes? let  lon2 = (lon+lon[i=@shf:1])/2
yes? let lat2 = if i ne 1 then (lat+lat[j=@shf:1])/2 ! Stop wrap.

yes? plot/dash/vs/ov/col=black/nolab lon2,lat2
yes? plot/dash/vs/ov/col=black/nolab transpose_xy(lon2),transpose_xy(lat2)

Andre, I think you want to be using the boundary specified version of the 3 argument shade plot if possible.


i.e. specify the cell boundaries for your shade plot command. You can then overlay those coordinates. That should fix the offset problem that you seem to be having.


On 10/09/15 10:29, Ansley Manke wrote:
If what you want is dashed lines, you could try fiddling with adding arguments to the /DASH qualifier, making shorter "down" segments:

Another alternative for showing the grids is to just make symbols, either at the grid cells or the corners,  so
yes? shade h,lon_hyc2,lat_hyc2
yes? plot/nolab/vs/ov/sym=dot lon_hyc2,lat_hyc2
or with a small symbol of some other shape.

That first plot you showed is pretty nice looking. Interestingly I can't seem to replicate what you did for a rectangular grid where I define 2-dimensional longitudes and latitudes for a simple 1-D grid, trying to simulate the 2-D coordinates:
yes? use coads_climatology

yes? let/title="longitude"/units=degrees_east lon = x[gx=sst] + 0*y[gy=sst]
yes? let/title="latitude"/units=degrees_north  lat = 0*x[gx=sst] + y[gy=sst]

yes? shade/l=1/hlim=180:240/vlim=0:40 sst, lon, lat
yes? let  lon2 = (lon+lon[i=@shf:1])/2
yes? let lat2 = (lat+lat[j=@shf:1])/2
yes? plot/dash/vs/ov lon2,lat2
I only see the horizontal dashed lines. Did you show your whole script?

On 9/9/2015 2:54 PM, Andre Paim wrote:
Dear all,

I've been trying to plot the grid from my model on top of a shade plot and I just cant seem to make it right. I have two grids, one regular and another irregular and curved along the coast. 

Here is how I'm plotting that for the regular grid.

   let lon_hyc2 = (lon_hyc+lon_hyc[i=@shf:1])/2
   let lat_hyc2 = (lat_hyc+lat_hyc[j=@shf:1])/2

   shade h,lon_hyc2,lat_hyc2
   plot/dash/vs/ov lon_hyc2,lat_hyc2

I don't know if that is a coincidence, but the plot works fine for the area I'm plotting as can be seem bellow. 
Inline images 3

However, when I try to do the same to the irregular grid, this is what I get.

Inline images 4

Both models are not projected, with lon and lat in their respective variables.
Is there an easier way to do that?

Thanks in advance,
Andre Rodrigues

[Thread Prev][Thread Next][Index]
Contact Us
Dept of Commerce / NOAA / OAR / PMEL / Ferret

Privacy Policy | Disclaimer | Accessibility Statement