[GRASS-dev] Extracting coordinates of RasterNumpy objects (pygrass)

[Pierre Roudier]

Sorry Pietro - more questions as I'm working my way through this new tool!

- As a follow-up from my last question, I am trying to convert the
pixels I have identified in my raster to a vector layer. While I
manage to create a Point() from scratch, I can't find out how I could
append coordinates to add more points to it. Using the previous
example:

>>> from grass import pygrass
>>> from pygrass.raster import RasterNumpy
>>> from pygrass.region import Region
>>> from pygrass.functions import pixel2coor
>>> elev = RasterNumpy("elevation")
>>> elev.open()
>>> a = elev > 144
>>> reg = Region()
>>> xcoords, ycoords = a.nonzero()
>>> # Switch to geographic space
>>>  coords = [pixel2coor(pixel, reg) for pixel in zip(list(xcoords), list(ycoords))]
>>>  # Add attribute column
>>>  xyz = np.column_stack((np.asarray(coords), z))
>>> from pygrass.vector import *
>>> p = pygrass.vector.geometry.Point() # Works
>>> p = pygrass.vector.geometry.Point(x = xyz[:,0], y = xyz[:,1], z=xyz[:,2]) # Fails

- Is there a method to create a set of points (class Point) and write
it into the GRASS database?

- I would like to use this {x,y,z} set of points in v.surf.bspline. I
don't suppose I got another choice but write my numpy array into the
GRASS database first right?

Thanks again for your work,

Pierre

2012/12/18 Pierre Roudier <pierre.roudier at gmail.com>:
> Thanks Pietro and Soeren,
>
> That's really great!
>
>
> I am looking forward to use this with v.surf.bspline to do on the fly
> interpolation (interpolation with SciPy is quite slow and memory
> hungry). I am coding an iterative process that calls the interpolation
> routine a lot of times (~10 to 40 interpolation steps), so I would
> like to make the most of GRASS interpolation power.
>
> Cheers,
>
> Pierre
>
> 2012/12/18 Pietro <peter.zamb at gmail.com>:
>> Hi Pierre,
>>
>> On Mon, Dec 17, 2012 at 9:44 PM, Pierre Roudier
>> <pierre.roudier at gmail.com> wrote:
>>> Thanks Pietro,
>>>
>>> Yes that answers the question, but just partly: I was actually
>>> wondering whether it would be possible to get the extracted x and y
>>> coordinates in geographic space (as opposed to the Numpy array space)?
>>
>> aah sorry, I forgot to say that you can use pixel2coord... see the
>> example below:
>>
>>>>> from grass import pygrass
>>>>> from pygrass.raster import RasterNumpy
>>>>> from pygrass.region import Region
>>>>> from pygrass.functions import pixel2coor
>>>>> elev = RasterNumpy("elevation")
>>>>> elev.open()
>>>>> a = elev > 144
>>>>> reg = Region()
>>>>> xcoords, ycoords = a.nonzero()
>>>>> for pixel in zip(list(xcoords), list(ycoords))[:5]:
>> ...     print pixel
>> ...
>> (0, 6)
>> (0, 7)
>> (0, 8)
>> (0, 9)
>> (0, 10)
>>>>> for pixel in zip(list(xcoords), list(ycoords))[:5]:
>> ...    print pixel2coor(pixel, reg)
>> ...
>> (228440.0, 630000.0)
>> (228430.0, 630000.0)
>> (228420.0, 630000.0)
>> (228410.0, 630000.0)
>> (228400.0, 630000.0)
>>
>> Best regards
>>
>> Pietro
>
>
>
> --
> Scientist
> Landcare Research, New Zealand



-- 
Scientist
Landcare Research, New Zealand