Opened 16 years ago

Closed 6 years ago

#2023 closed enhancement (fixed)

Add a Python wrapper for GDALGridCreate() function

Reported by: dron Owned by: hobu
Priority: normal Milestone:
Component: PythonBindings Version: unspecified
Severity: normal Keywords:
Cc: warmerdam, chrisbarker, Ari Jolma

Description

The recently added GDALGridCreate(http://www.gdal.org/gdal__alg_8h.html#272e099fee92569aff154bbb035fc7c1) function needs to be wrapped for use in Python. The function takes three floating point arrays pdfX, pdfY and pdfZ as input and outputs an interpolated data as a regular array of specified data type and dimensions.

To control the interpolation process function takes the algorithm name and algorithm options. There are different set of options for every particular algorithm. At the moment we have GDALGridInverseDistanceToAPowerOptions (http://www.gdal.org/structGDALGridInverseDistanceToAPowerOptions.html), GDALGridMovingAverageOptions(http://www.gdal.org/structGDALGridMovingAverageOptions.html) and GDALGridNearestNeighborOptions(http://www.gdal.org/structGDALGridNearestNeighborOptions.html) structures for three available algorithms. I am thinking about passing this options to Python function in form of dictionary:

options["algorithm"] = "nearest"
options["radius1"] = 12.3
options["radius2"] = 34.5

The whole Python API could look like that:

GridCreate(numpy_array, (minx,miny,maxx,maxy), nXSize, nYSize, eType, options)

numpy_array is a 2D NumPy array containing all three pdfX, pdfY and pdfZ vectors.

Change History (9)

comment:1 by hobu, 16 years ago

Cc: warmerdam added

CC'ing Frank for any input he might have on this

Here was our IRC output:

DronK>	I would like to wrap this function: http://www.gdal.org/gdal__alg_8h.html#272e099fee92569aff154bbb035fc7c1
	<DronK>	But I am not confident with the new style bindings.
	<hobu>	File a bug for me and I'll do it
	<hobu>	what is void* pData?
	<DronK>	That would be cool!
	<hobu>	inputted data, or outputted?
	<DronK>	Outputted.
	<hobu>	I am working on the progress callbacks right now, and I'm very close to having them work right
	<DronK>	This is array, controlled by the nXSize/nYSize/eType
	<hobu>	This will get us back to where we were with the old bindings with respect to the gdal operations
	<hobu>	array == Band? or Dataset?
	<even_>	DronK: if it's new API, poOptions, pdfX, Y and Z could be "constified"
	<DronK>	Well, array is an array here. This particular function does not recognize the band or dataset concepts.
	<DronK>	even_, good point.
	<even_>	but we won't start the const battle as Mateusz would like ;-)
	<DronK>	Happy me, I have missed that battle. But I am totally agree that these ones should be cons pointers.
	<even_>	If I remember well, Frank was a bit worried that Mateusz would start a const battle in GDAL
	<hobu>	The list of points should be just a dumb list of floats?
	<hobu>	also, will the function default its window to the extent of the points if the dfXMin, etc are 0?
	<hobu>	also, what is poOptions, and how do you foresee creating that from Python?
	<DronK>	hobu, yes, the function gets three floating point arrays (that should be Numeric arrays I think) and returns back regular grid.
	<DronK>	dfXMin and other must be set, there is no default, it is essential parameter.
	<DronK>	pOption is a tricky thing here. I will be more verbosive in ticket.
	<hobu>	I think they should be tuples or lists --- [[1,2,3],[4,5,6],[7,8,9]]
	<hobu>	We'd have to do something really special to get numpy arrays
	<DronK>	They can be large. Really large. Will be lists effective for that?
	<hobu>	Not really :)
	<hobu>	Let me talk to Chris Barker though
	<DronK>	Conversion fom list to C array could be not effective also.
	<hobu>	IIRC, numpy has the ability to give you a handle to its internal array without copying
	<hobu>	But I don't know the particulars
	<DronK>	Yes, that is why I want numpy here :-)
	<hobu>	Additionally, we can't get too fancy, because we have no numpy included in the general bindings (it is an additional, stand alone module)
	<DronK>	Is it possible to build this binding when numpy is switched on only?
	<hobu>	I think we could have the wrapper be built by default, and then use the swig typemap system to do the argument coercion magic
	<hobu>	this way the other languages can come up with their strategies for dealing with this function
	<hobu>	An issue is the typemaps stuff does not know about numpy
	<hobu>	knowing about numpy == including numpy's headers.... this can be tough to do
	<hobu>	I could try to fetch the numpy headers from ./configure, however, but we've never gone that far.
	<hobu>	Also, what to do when we can't include numpy headers
	<hobu>	(we've always supported building the bindings without numpy/Numeric)
	<DronK>	Will it be possible to extend this API in the future? We can go with lists now and add numpy support later if it is poddible.
	<DronK>	The main reason for me to get this binding working now is a test suite.
	<hobu>	I would propose that we make numpy *required* to build the next-gen bindings. But this might take a PSC motion
	<DronK>	I am on that, certainly.
	<hobu>	Additionally, GridCreate is a very generic function name... a sufficiently brainwashed ESRI user might think it is the function to open up a dataset :)
	<hobu>	API like so? GridCreate(numpy_array, (minx,miny,maxx,maxy), nXSize, nYSize, eType, options)
	<DronK>	Huh, I even wrote some docs, so...
	<hobu>	Python users don't read docs :P
	-->|	epifanio (n=chatzill@host21-80-dynamic.53-82-r.retail.telecomitalia.it) has joined #gdal
	<DronK>	I think that we can pass three 1D arrays instead of one 2D, it looks simpler for me.
	<DronK>	Options should be a dictionary.
	<hobu>	I was thinking numpy_array was a 3D array, and the user would be responsible for organizing it how they wished
	<hobu>	but I'm flexible either way
	<hobu>	poOptions is a CSLString?
	<DronK>	I think is is depends on taste. I like three arrays, but I have nothing against single array. Note, that with three arrays we can have 4GiB data in them.
	|<--	epifanio has left epimetheus.hobu.net:57000 (Client Quit)
	<hobu>	If the users passes in an array that isn't nPoints long, you'll have trouble iterating
	<hobu>	if it is a single array, we can do a check in the typemap to ensure that things are in alignment
	<DronK>	Hmm, that right.
	<DronK>	Let it be a single array then.
	<CIA-17>	rouault * r13021 /trunk/gdal/ogr/ogrsf_frmts/gpx/ogrgpxlayer.cpp: Support properly <time> content as being OFTDateTime
	<hobu>	If the user needs to do *big* data, running this little wrapper in Python isn't really for them :)
	<CIA-17>	rouault * r13022 /trunk/autotest/ogr/ (data/test.gpx ogr_gpx.py): Update GPX test to test dateTime
	<DronK>	poOptions is not a string, it is an ugly thing. It is a pointer to an object of several possible types: http://www.gdal.org/structGDALGridInverseDistanceToAPowerOptions.html
	<DronK>	Or this one: http://www.gdal.org/structGDALGridMovingAverageOptions.html
	<nhv>	DronK are you familiar with this work http://www.marine.csiro.au/~sak007/
	<hobu>	yikes
	<hobu>	what about krig and spline? Would they have their own options structs?
	<DronK>	What is exactly used depends on the first parameter.
	* DronK	thinks he scares hobu.
	<hobu>	I'm not scared of wrapping it (it looks fairly straightforward except for the numpy array stuff).
	<hobu>	I'm scared of how is a user going to use it effectively
	<DronK>	nhv, actually I didn't used that code, and I am only implemented simple things for now.
	<hobu>	(anything we put in the bindings, people end up using, and this causes us misery :) )
	<nhv>	cool
	<DronK>	nhv, I am planning to add more algorithms, but I need more general purpose math to be implemented in GDAL kernel first.
	<nhv>	gottcha
	<hobu>	our discussion should have enough info for a ticket. I will try to get to it after I finish callbacks today
	<DronK>	hobu, that was the simplest way to do in C. Every algorithm has its own set of options and we can constanly increase the number of supported algorithms not changing the C API.
	<hobu>	yes, I understand. It's a complex operation. I wasn't saying you made it overly complicated :)
	* nhv	wonders how much 'grid math' gdal really wants to have, historically it has been about data access not image processing
	<DronK>	I thought about dictionary (something like that in OGRStyle API), but this requres adding option init/delete routines... So I have drop that.
	<hobu>	I expect I will make some class(es) to behave correctly in Python for the options
	<DronK>	nhv, no rocket scuience math, just least squares/fitting routines, we already have these ones.
	<DronK>	The algorithm name could be added to options in Python, not the separate parameter like in C.
	<nhv>	understood but once started on a path it seems the road turns into a highway
	<DronK>	He-he, no vehicle to ride that higway.

comment:2 by hobu, 16 years ago

Andrey,

I don't think I will be able to get to this before the release on Tuesday...

I think we would want to use the numpy swig typemap stuff from numpy's subversion repository to make the array coercion easier http://svn.scipy.org/svn/numpy/trunk/numpy/doc/swig/

comment:3 by hobu, 16 years ago

Cc: chrisbarker added

cc'ing Chris in case he can add some insight to what's going on in swig/numpy land, as I've seen some traffic on the swig lists lately...

comment:4 by dron, 16 years ago

Howard,

There is no any urgency with this wrapper. I just need it to simplify a tests creation before I will continue GDALGrid development and that will go in new post-1.5 GDAL branch anyway.

Best regards,
Andrey

comment:5 by Ari Jolma, 16 years ago

Cc: Ari Jolma added

Hmm, this would of course be useful for Perl wrappers too. But.. in fact I'd like the API to be something like dataset = GridCreate(layer_of_multiple_points, ...), so it would be easy to integrate into an application.

But nevermind, I can build the extra layer in Perl bindings too.

Regards,

Ari

ps: isn't this an enhancement and not a defect?

comment:6 by dron, 16 years ago

Type: defectenhancement

Ari,

Do you mean to pass OGR layer as an input parameter? I would like not to do that for four reasons:

  1. It is better to keep bindings for different languages with the same set of parameters.
  1. Often you need to process the SQL query results and not a single layer, so functionality of the GDALGridCreate() should be flexible enough to allow this.
  1. Sometimes you want to get Z value from the semantic field and not from the geometry, so we need more flexible function again.
  1. It is very common to pass XYZ data in form of simple CSV file, that is very easy to parse in script and create input array without any connection to OGR at all.

If you really like a function that gets OGR layer as input I can create an appropriate C wrapper on top of GDALCreateGrid() and that function could be exposed in bindings along with GDALCreateGrid().

Best regards,
Andrey

comment:7 by wyogeo, 10 years ago

I just ran across this because I find myself in need of this functionality. Was this enhancement ever done? I've searched through the gdal import and can't find a GridCreate function.

comment:8 by Even Rouault, 10 years ago

The Java bindings has been done, but the Python bindings still not. Contributions or funding are always welcome.

comment:9 by Even Rouault, 6 years ago

Resolution: fixed
Status: newclosed

gdal.Grid() is available in Python bindings

Note: See TracTickets for help on using tickets.