[gdal-dev] How to use TransformPoints in python?

[Hogendoorn, Rene]

The following simple python script

import os, osgeo, numpy
from osgeo import ogr, osr


def getProjection( EPSGcode ):
    """ set a projection from ... to lat-long """

    srs = osgeo.osr.SpatialReference()
    srs.ImportFromEPSG( 4326 ) # WGS84
    ll = srs.CloneGeogCS()

    srs.ImportFromEPSG( EPSGcode )
    # RD-NL lacks TOWGS84
    if EPSGcode == 28992:
        srs.SetTOWGS84( 565.237, 50.0087, 465.658, -0.406857, 0.350733, -1.87035, 4.0812 )

    return osr.CoordinateTransformation( srs, ll )


if __name__ == '__main__':

    proj = getProjection( 28992 )

    polygon =[[132817.006604435708141, 550302.852720651309937, 0.],
              [131182.28895997320069, 558340.214472591876984, 0.],
              [132578.61028128489852, 558748.893883707583882, 0.],
              [136631.347774848196423, 553436.061539204441942, 0.],
              [136631.347774848196423, 553436.061539204441942, 0.],
              [132817.006604435708141, 550302.852720651309937, 0.]]

    x = numpy.array( [p[0] for p in polygon] )
    y = numpy.array( [p[1] for p in polygon] )
    z = numpy.array( [p[2] for p in polygon] )

    print proj.TransformPoints( len(polygon), x, y, z)

returns:

Traceback (most recent call last):
  File "<stdin>", line 37, in <module>
  File "/var/tmp/libgdal0-1.6.0-build/usr/lib64/python2.6/site-packages/osgeo/osr.py", line 622, in TransformPoints
    
TypeError: in method 'CoordinateTransformation_TransformPoints', argument 3 of type 'double *'

I also tried straight python lists of doubles.
So what type of arguments does TransformPoints expect (as there is obviously nothing like double* in python)?

Regards

Rene
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://lists.osgeo.org/pipermail/gdal-dev/attachments/20090116/da2f8cec/attachment.html