obs.x (and all others) wtr to python start 0 and fortan (ROMS) start 1?

[ivicajan]

Brian & Dale, I am puzzled as when you call seapy gridder and use grid.ij (or grid.ijk) you have all indexes starting from 0 (as python counts)/ You need that when searching in the matrices etc. However, when putting all into netcdf we have to follow fortran which use start from 1. Is it possible that seapy created obs.x, .y, .z which are -1 from the one created before with matlab's obs_gridder? Should we in to_netcdf add obs.x += 1 (the same for y)? Those are the only indexes, others are real locations, values...

Cheers, Ivica

[powellb]

It would appear that we are missing a (+1) to the ii and jj fields in the gridder before the obs object is constructed.

grid.ij works as expected:

j,i = grid.ij(grid.lon_rho[10,10], grid.lat_rho[10,10]))

should return 10,10. That is working; however, as you note, the obs should count from 1, so that is a bug in the gridder. Will investigate a bit more and upload change.

[ivicajan]

Agree that gridder is working fine, but in the last step when you write to netCDF should follow ROMS convention and add +1 for index variables only. If you mess with obs structure then within python you will not be consistent.

[powellb]

The obs structure is always a ROMS construct and should always count from 1, not 0. to_netcdf is not where a +1 would happen because everytime you called it the obs would shift.

I will make the changes tomorrow because there is an impact elsewhere.

On Sep 7, 2016, at 9:39 PM, ivicajan notifications@github.com wrote:

Agree that gridder is working fine, but in the last step when you write to netCDF should follow ROMS convention and add +1 for index variables only. If you mess with obs structure then within python you will not be consistent.

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

[ivicajan]

For example, when I want to kill all obs in the shallower then say 10m then I take info from obs.x and obs.y and they should be consistent with python way of arrays: x = np.round(obs.x).astype(int) y = np.round(obs.y).astype(int) h = grid.h[y,x] # my grid is having start at 0 index idx = np.where(h < 10)

[powellb]

Hmmmm, I need to go back through the original procedure. The python behaves exactly as the previous matlab routines did (see 'latlon_obs.m'). We would be having masking issues if there was a +1 offset, which we don't see. I need to go through old notes tomorrow, but for now, it seems correct.

On Sep 7, 2016, at 9:50 PM, ivicajan notifications@github.com wrote:

For example, when I want to kill all obs in the shallower then say 10m then I take info from obs.x and obs.y and they should be consistent with python way of arrays: x = np.round(obs.x).astype(int) y = np.round(obs.y).astype(int) h = grid.h[y,x] # my grid is having start at 0 index idx = np.where(h < 10)

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

[powellb]

I am going to close this and say it is correct as is.

Refer to the ROMS source extract_obs.F. In the top description, it states:

! All the observations are assumed to in fractional coordinates with ! ! respect to RHO-points: !

there is shown a map of points, and the rho-points start counting from 0 (with u and v starting from 0.5). This is how the original matlab-based observational routines were done (starting from 0) as well.

On Sep 7, 2016, at 9:50 PM, ivicajan notifications@github.com wrote:

For example, when I want to kill all obs in the shallower then say 10m then I take info from obs.x and obs.y and they should be consistent with python way of arrays: x = np.round(obs.x).astype(int) y = np.round(obs.y).astype(int) h = grid.h[y,x] # my grid is having start at 0 index idx = np.where(h < 10)

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

[ivicajan]

No problem, I was comparing something previously done in matlab and now in seapy and found mismatch, basically using index in matlab and those in the python. The origin of all problems was start counting of arrays from 0 or 1.

If you are saying that ROMS and python are counting fractional index from 0 then if I used obs_read.m (that is still 0 based) and then I load gridwith grid_read.m and have for example depths grid.h(round(obs.y), round(obs.x)) it will be off to lower-left corner (in the matlab as matlab matix is from 0)?

Then obs.z are exception and they should go from 1 to N? Remember we used before surface and put number of levels to obs.z (for say SST)..

Cheers I

On 08/09/16 16:07, Brian Powell wrote:

I am going to close this and say it is correct as is.

Refer to the ROMS source extract_obs.F. In the top description, it states:

! All the observations are assumed to in fractional coordinates with ! ! respect to RHO-points: !

there is shown a map of points, and the rho-points start counting from 0 (with u and v starting from 0.5). This is how the original matlab-based observational routines were done (starting from 0) as well.

On Sep 7, 2016, at 9:50 PM, ivicajan notifications@github.com wrote:

For example, when I want to kill all obs in the shallower then say 10m then I take info from obs.x and obs.y and they should be consistent with python way of arrays: x = np.round(obs.x).astype(int) y = np.round(obs.y).astype(int) h = grid.h[y,x] # my grid is having start at 0 index idx = np.where(h < 10)

— You are receiving this because you commented. Reply to this email directly, view it on GitHub, or mute the thread.

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/powellb/seapy/issues/29#issuecomment-245524137, or mute the thread https://github.com/notifications/unsubscribe-auth/ARvJsGkeS_ShEvMtG2s1la3CLntiHgkwks5qn8IngaJpZM4J3o2k.