[GRASSLIST:1083] Reproject raster data - weighted average

Hi,

I have some floating point raster data sets in a lat/long location. I want to reproject these data into a cylindrical equal area location. Obviously I can do this using `r.proj' but I want the values of cells in the reprojected map to be calculated as averages of the corresponding cells in the original projection, weighted by their relative area in the new grid cells.

Is there any way to do this in Grass?

Thanks,
David

Hello
Not totally clear on what you mean if 'r.proj method=bilinear' does not work
for you, what about convering the original raster to sites with r.to.sites
and then re-projection with s.proj, and use one of the various
interpolation methods (e.g. s.surf.idw or s.surf.rst) to get the raster
you want in the new location? Appropriate solution sort of depends on the
relative cell resolutions of the two locations, which you didn't mention.

Paul

On Tue, 2 Sep 2003, David Orme wrote:

Hi,

I have some floating point raster data sets in a lat/long location. I
want to reproject these data into a cylindrical equal area location.
Obviously I can do this using `r.proj' but I want the values of cells
in the reprojected map to be calculated as averages of the
corresponding cells in the original projection, weighted by their
relative area in the new grid cells.

Is there any way to do this in Grass?

Thanks,
David

David,

You may try to use r.statistics to compute average values before you
reproject your data (method=average) Something like:

1. Prepare new map (X) in your lat./long. location with the desired
resolution, lower than original - larger pixel size (g.region, r.mapcalc
may be used)

2.Run in your lat./long.
r.statistics cover=Y base=X method=average output=N

N - new map with average of Y values based on X.

3. r.proj

Cheers,

Rado

Dòa Ut, 2003-09-02 at 11:08, David Orme napísal:

Hi,

I have some floating point raster data sets in a lat/long location. I
want to reproject these data into a cylindrical equal area location.
Obviously I can do this using `r.proj' but I want the values of cells
in the reprojected map to be calculated as averages of the
corresponding cells in the original projection, weighted by their
relative area in the new grid cells.

Is there any way to do this in Grass?

Thanks,
David

--
Radoslav Bonk M.S.
Dept. of Physical Geography and Geoecology
Faculty of Sciences, Comenius University
Mlynska Dolina 842 15, Bratislava, SLOVAKIA
tel: +421 905 968 127 e-mail: rbonk@host.sk

Hi,

Providing more information... I need to transform some global data sets in geographic resolutions, for example quarter degree, into a cylindrical equal area grid. The resolution of the CEA location has been chosen so that it matches the geographic resolution of the data set at the latitude of true scale - as a result the longitudinal divisions of the two grids coincide but the latitudinal lines divisions differ (except at the equator - the origin of both grids). What I would like to be able to do is to calculate the values of the new grid as an average of the original cells underlying each new cell weighted by their contribution to the area of the new cell. As a toy example, say I have three cells going north from the equator with values 5, 7, 8 and the new grid happens to form two cells with a division half way through the second cell. The first cell in the new grid should contain the mean of 5 and 7 weighted by 1 and 0.5 (i.e. 5.6667) - the second the weighted mean of 7 and 8 weighted by 0.5 and 1 (7.6667).

Paul and Rado (thanks for the amazingly fast replies!) suggested:

1) Interpolation from a point surface - I think this would give the correct answer if the interpolation neighbourhood was set up correctly. Unfortunately because I'm going to an equal area grid, the neighbourhood for cells in the new grid varies with latitude.

2) r.neighbours - again I think this would give the correct answer locally but because the CEA and geographic grids do not overlie each other in a regular fashion it won't work beyond that, unless there is a way of feeding a grid with varying North South resolution into r.neighbours, which I seriously doubt.

I can think of a way of doing it using a transformation matrix in R, for instance, but I wondered if I was reinventing the wheel.

Thanks,
David

On Wed, 3 Sep 2003, David Orme wrote:

Hi,

Providing more information... I need to transform some global data sets
in geographic resolutions, for example quarter degree, into a
cylindrical equal area grid. The resolution of the CEA location has
been chosen so that it matches the geographic resolution of the data
set at the latitude of true scale - as a result the longitudinal
divisions of the two grids coincide but the latitudinal lines divisions
differ (except at the equator - the origin of both grids). What I would
like to be able to do is to calculate the values of the new grid as an
average of the original cells underlying each new cell weighted by
their contribution to the area of the new cell. As a toy example, say I
have three cells going north from the equator with values 5, 7, 8 and
the new grid happens to form two cells with a division half way through
the second cell. The first cell in the new grid should contain the mean
of 5 and 7 weighted by 1 and 0.5 (i.e. 5.6667) - the second the
weighted mean of 7 and 8 weighted by 0.5 and 1 (7.6667).

Paul and Rado (thanks for the amazingly fast replies!) suggested:

1) Interpolation from a point surface - I think this would give the
correct answer if the interpolation neighbourhood was set up correctly.
Unfortunately because I'm going to an equal area grid, the
neighbourhood for cells in the new grid varies with latitude.

2) r.neighbours - again I think this would give the correct answer
locally but because the CEA and geographic grids do not overlie each
other in a regular fashion it won't work beyond that, unless there is a
way of feeding a grid with varying North South resolution into
r.neighbours, which I seriously doubt.

I can think of a way of doing it using a transformation matrix in R,
for instance, but I wondered if I was reinventing the wheel.

Were it not for your wanting a weighted average, I think gdalwarp would do
this, but it is restricted to nearest neighbour resampling, as I
understand.

Roger

Thanks,
David

--
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Breiviksveien 40, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 93 93
e-mail: Roger.Bivand@nhh.no