Hi Chris and Jaison,
I have just one thing to add to what looks like a great answer from
Jaison.
Chris mentioned in his original question that the SCAT2GRID function
smears
out the results onto land. You can try choosing the SCAT2GRID
parameters
to get a closer fit, or use a mask on the final result.
If var_2 is the data on the X-Y-T grid, then you could do something
like this:
... compute var_2
use etopo20 ! choose a data set with a resolution similar to your
output grid
! Regrid the elevations to your grid and use them to define a mask
which is
! 1 over the ocean and missing everywhere else
let land_mask = if rose[d=2] LT 0 then 1
let masked_var = land_mask * var_2[d=1]
shade masked_var
Ansley Manke
Jaison Kurian wrote:
Hi Chris,
You reached there almost...but missed few points. I modified
your script (see the attachment:scat2grid.jnl, using the data
sample_data.txt) to get proper time dimension for the variable based on
following details...
1. Your data is repeating for the common X-Y locations for 5 time
points ( 1990- 3/31, 4/30, 5/21, 6/14 & 6/21). But the lat-lon
point "80.50 -67.75" exists only for "1990 03 31". If this
point is skipped we are left with 107 spacial points and 5 time
points ==> 107*5 = 535 values (or lines in the data file).
Now read the variable in Ferret like...
define axis/x=1:107:1 xscat ! spatial points
define axis/t=1:5:1 tscat ! time points
define grid/x=xscat/t=tscat gscat
file/grid=gscat/var=..............
2. to define a time axis choose only the five values from year, month
& day (skip all the repeated values). Get the time values using
days1900 function, define a time axis..assign to the variable
using @ASN regridding.
3. Use this new variable in scat2grid function. Now the variable is
2D but one of the dimention is time. So you have to give only
the common lat and lon points as the arguments to the scat2grid.
Define the final destination grid such that it is within the limits
of the original data. Dont forget to use the /UNITS attribute
while defining latitude or longitude axis
4. Define a climatological monthly time axis then regrid to this
usin @MODULO regridding to get the climatology.
Please go through the modified script..if anything need to be
explained..let me know.
As a "bonus" you can see how to use "POLYGON" command to look at the
original data in a spatial plot without any interpolation/extrapolation.
Hope this helps
With Regards
Jaison
On Thu, 18 Nov 2004, Chis Konig wrote:
Hi there,
I have data that came in 2D-fields at irregular times, i.e. it is
unevenly spaced in time (weekly - monthly) and space (i.e. 2D but was
not on a lon-lat-grid). I'd like to make a climatology with monthly
means etc. from it.
Unless there is an elegant way of doing that all at once (3D
regridding?) I think I have to achieve this in two steps (but correct me
if I'm wrong):
- First regrid it in space onto a lon-lat-grid for each time step.
- Secondly regrid it over time, i.e. make monthly averages etc.
I might prefer a two-stepped approach anyway as it gives me more
control, e.g. different types of interpolation for space and time.
The data is in ASCII and in the form:
yr mo dy lat lon var
An example file can be found here:
http://math.nyu.edu/~konigc/sample_data/sample_data.txt
The first three variables (time) are non-decreasing but identical for
several hundreds to thousands data points (all the data from one field).
That's how far I got:
!!!!!!!!!!!!!!!!!!!!!!!!
define axis/x=1:1000:1 x_temp
define grid/x=x_temp g_temp
file/var="yr,mo,dy,lat,lon,var"/grid=g_temp sample_data.txt
let gregday = days1900(yr,mo,dy)
define axis/x=-75:-55:0.25 xax
define axis/y=79:83:0.25 yax
!Those axes would be different for the total data.
let var_2=scat2gridgauss_XY(lon,lat,var,x[gx=xax],y[gy=yax],0.5,0.5,2,2)
!OR
let var_2=scat2gridlaplace_XY(lon,lat,var,x[gx=xax],y[gy=yax],5,2)
shade var_2
go land
!!!!!!!!!!!!!!!!!!!!!!!!
My problems with that approach:
- I lose the time dimension in the regridding process.
- The data is confined to the ocean. But with scat2gridXXX the
information 'spills' onto the land. Is there a possibility to use a
simpler interpolation (nearest neighbour, linear...)?
Thanks for any help!
Chris
!
! delete the line containing data for 1990 03 31 80.50 -67.75
! to bring all days to a common scattered points
!
define axis/x=1:107:1 xscat ! data points
define axis/t=1:5:1 tscat ! time points
define grid/x=xscat/t=tscat gscat
file/var="yr,mo,dy,lat,lon,var"/grid=gscat sample_data.txt
! let us have a cross check
list/y=5 var
! get all the 5 dates
let yr_5 = yr[i=1] ! values need only for I=1,L=1:5
let mo_5 = mo[i=1] ! "
let dy_5 = dy[i=1] ! "
list yr_5, mo_5, dy_5
let gregday = days1900(yr_5,mo_5,dy_5)
define axis/t/t0="01-JAN-1900:00"/units=days/from_data t5day=gregday
list t[gt=t5day]
! get the variable onto the time grid instead of Y
let var_t = var[GT=t5day@ASN]
plot/l=5 var_t
message
! now do the scat2grid
! choose min & max for defining the axes so that those
! values are within the actual range
let lat_1D = lat[l=1] ! we need only the 107 common values
let lon_1D = lon[l=1] ! we need only the 107 common values
define axis/x=-61.5:-67.5:0.25/units=longitudes xlon
define axis/y=80.5:81.75:0.25/units=latitudes ylat
let var_2=scat2gridgauss_XY(lon_1D,lat_1D,var_t,x[gx=xlon],y[gy=ylat],0.125,0.125,2,2)
fill/l=1 var_2
plot /ov/vs/color=7 lon[l=1],lat[l=1]
go land
! no you can define a monthly climatological axis (if needed..because now
! you have one value per month)..and make the climatology using var_t
! let us see what "polygon" can do
message
let xpts = XSEQUENCE(lon[l=1,i=1:107])
let ypts = XSEQUENCE(lat[l=1,i=1:107])
let dxp=0.2 ; let bx=ysequence({`(-1)*dxp`,`dxp`,`dxp`,`(-1)*dxp`})
let dzp=0.05 ; let bz=ysequence({`dzp`,`dzp`,`(-1)*dzp`,`(-1)*dzp`})
polygon bx+xpts,bz+ypts,var_t[l=1]
plot /ov/vs lon[l=1],lat[l=1]
go land
|