How can I get the pixel data corresponding to a shapefile with some coordinates' points in IDL?

Question

What I want to get is the pixel's data corresponding to these points in the vector file. And add the maximum values arround the pixel in a window of 8 to the coordinate as a new value. I am a newcomer to IDL, how can I solve the problems ? Thank you in advance!

Answer 1

I'm not quite sure what pixel data you're referring to, so I'll assume it's another array which you can index with the shapefile coordinates. Perhaps this article will help you to extract the vertices of the polygon (looks like you'll have to cough up a bit of cash for the routine he uses to do so). If you want to get the maximum values around each pixel, after you get each x and y vector to index into the data (pixel below) you'll have to do some dimensional juggling. Something like the following might work, but I haven't tested it (you'll probably have to play with the dimension over which the maximum is calculated):

n = n_elements(x)
; pixel = findgen(np, np)

; Make a copy of "pixel" so we can eliminate the center element
; from consideration of the maximum
new_pixel = pixel
new_pixel[x, y] = new_pixel[x + 1, y + 1]

x = rebin(reform(x, n, 1), n, 9)
y = rebin(reform(y, n, 1), n, 9)

x += rebin(reform([-1, 0, 1], 1, 3), n, 9)
y += rebin(reform(rebin(reform([-1, 0, 1], 3, 1), 3, 3), 1, 9), n, 9)

; Get surrounding elements
new_pixel = new_pixel[x, y]
; Get maximum
max_pixel = max(new_pixel, dimension=2)

Answer 2

I highly recommend using QGIS for shapefile operations before importing the matrix indices into IDL. Then you can perform a simple calculation:

output = input[where(index eq shp_value)]

To align the pixels, I would write the shapefile values to an existing raster when creating the index raster in QGIS. I believe you want the Rasterize tool. For more information, read my reply regarding zonal statistics.

Cheers,

Adam