stations output for u and v seems to be shifted

[wmartin]

I seems to me that the interpolation code for stations is shifting u and v by one node. For example, if you have a uniform x grid at 1km spacing, the value of u at the 8km station seems to be the average of values at 6.5 km and 7.5 km in the u grid (rather than 7.5 and 8.5). I see a similar thing for v in the y direction. I have STATIONS_CGRID undefined and all my stations are specified with interger grid pairs like 7,8 with FLAG = 0. My grid has both positive and negative values of x and y but the problem occurs in all quadrants. I found this by comparing data directly from the history files with that in the station files and I have checked the values of x_rho, x_v, x_u etc. from those files to try and make sure I'm not just shifting something myself. For example the values from history for u come from the grid points that have x_u = 6.5 and 7.5 which are where I expect them. I haven't dived into the code, hoping that someone can either explain what I'm doing wrong or has already got the fix.
Wayne

[klyons]

I have a question on a similar vein regarding the provision of X-POS and Y-POS in a stations.in file in terms of I and J grid pairs. Am I correct in thinking that the I and J pairs provided in the stations.in file are used as zero-based indices for data arrays in ROMS (e.g. grid size in Rho points is 0:Lm+1 by 0:Mm+1)? Also, are the exterior walls (0, Lm+1 ,0, Mm+1) considered invalid for station positions?

Does this also apply to the floats.in file?

Does it apply to position of point sources?

Thanks,
Kieran

[msd]

I just discovered the same problem with the station output for u and v being shifted (w/ STATIONS_CGRID undefined), but I think the problem is only for the 3d fields and the fix is pretty simple. Take a look at where interpolation is done for u-points in extract_sta.F (version 2.2, don't know if this has been fixed in 3.0).

For the 2D fields:

Code: Select all

      ELSE IF (gtype.eq.u2dvar) THEN
        DO np=1,Npos
          Xgrd=Xpos(np)+Xoff
          Ygrd=Ypos(np)
          i1=INT(Xgrd)
          i2=i1+1
          j1=INT(Ygrd)
          j2=j1+1

so given Xpos(np) = 1, i1 = 1 and i2 = 2 which I think is OK for the u-grid indices. However, for the 3D fields:

Code: Select all

      ELSE IF (gtype.eq.u3dvar) THEN
        DO np=1,Npos
          Xgrd=Xpos(np)-Xoff
          Ygrd=Ypos(np)
          Zgrd=Zpos(np)
          i1=INT(Xgrd)
          i2=i1+1
          j1=INT(Ygrd)
          j2=j1+1

which, given Xpos(np) = 1, gives i1 = 0 and i2 = 1 which is not OK for the u-grid indices bracketing a rho-grid point at i=1 (and would shift things in the direction Wayne found).

If this is true, then the fix should just be to change '-Xoff' and '-Yoff' to '+Xoff' and '+Yoff' in extract_sta3d in the extract_sta.F file.
Mike