Kenneth Bowman wrote:
> In article <1102591840.222724.104440@f14g2000cwb.googlegroups.com>,
> "Martin" <m.doyle@uea.ac.uk> wrote:
>
>
>> Hello everyone,
>>
>> I've been trying to find a way of doing this for ages, but my head is
>> sore after banging it against my computer screen!
>>
>> I have an irregularly gridded dataset (gaussian grid) and I want to be
>> able to 1) interpolate it onto a regular grid and 2) plot the data so
>> that each grid square is coloured depending on the value within each
>> square.
>>
>> I have a problem with 2) in that I can't find any IDL routines that do
>> this. Does anyone happen to know how to go about this at all? I'd be
>> grateful for a any suggestions you might have.
>
>
> There are two ways to do this that come to mind.
>
> One is to use POLYFILL to draw a polygon for each grid cell with the
> desired color. This is probably not the best way, as the PS device only
> supports 256 colors for line graphics, and it may look ugly on some map
> projections (e.g., when grid cells are not rectangles).
>
> The other way is to create an image, but not interpolate the data:
>
> MAP_SET, /HAMMER, LIMIT=limit, /ISOTROPIC, /NOBORDER
>
> data = DIST(15)
> bilinear = 0
> proj = MAP_IMAGE(data, i0, j0, ni, nj, $
> COMPRESS=1, BILINEAR=bilinear, $
> LATMIN=-90.0, LONMIN=0.0, LATMAX=90.0, LONMAX=360.0)
>
> image = BYTSCL(proj)
> IF (!D.NAME EQ 'PS') THEN $
> TV, image, i0, j0, XSIZE=ni, YSIZE=nj $
> ELSE $
> TV, image, i0, j0
>
> MAP_CONTINENTS
> MAP_GRID, GLINESTYLE=0
>
>
> Change bilinear to 1 to see the effects of interpolating. You can do
> something more complicated than a simple BYTSCL to define the colors for
> your data.
>
> Ken Bowman
Or, if the data sits on a nice rectangular grid, you could just tv a rebinned
version of the array over a plot which establishes your coordinates
system, like this:
;input data
im=dist(10,10)
xrange=[-2.,2]
yrange=[-5.,5]
;establish coordinates
plot,[0,0],[0,0],/nodata,/xstyle,/ystyle,xrange=xrange,yrang e=yrange
;rebin the data array to the right size
res=convert_coord(xrange,yrange,/data,/to_device)
xsize=round(res[0,1]-res[0,0])
ysize=round(res[1,1]-res[1,0])
im2=congrid(im,xsize,ysize,/center)
;plot the image
tv,bytscl(im2),-2,-5,/data,xsize=xsize,ysize=ysize
;overplot of the axis (since the tickmarks were covered by the image)
plot,[0,0],[0,0],/nodata,/xstyle,/ystyle,xrange=xrange,yrang e=yrange,/noerase
This is fast, but works just for device with pixels...
Cheers,
Paolo
--
____________________________________________________________ ________
Paolo Grigis
ETHZ - Institute of Astronomy email: pgrigis@astro.phys.ethz.ch
Scheuchzerstrasse 7
ETH Zentrum phone: ++41 1 632 42 20
8092 Zurich fax : ++41 1 632 12 05
Switzerland http://www.astro.phys.ethz.ch/
____________________________________________________________ ________
|