[gdal-dev] CS transformations

[Al Niessner]

I am a newbie at osgeo and gdal and have done the obligatory googling
and reading of tutorials and documentation. I have made it pretty far
before I found my way here so thanks to everyone for all of the good
documentation and references online already.

language: Python 3.4
OS: Ubuntu 14.04 (up to date)

Task:
Read in various georef PDF files, make a bounding box of the valid data
in the pixel or image coordinate system (CS) and WGS84. The data
somewhat mimics the USGS in that there is a metadata item called
NEATLINE that describes the valid data polygon. So, I thought I would
use a USGS historical quad to build up my software decoding of the maps.
I made it this far...

# read the PDF
ds = gdal.Open('Madulce_Peak_O34119F5_geo.pdf')

# get the affine to transform to pixel CS
at = ds.GetGeoTransform()

# get the dataset projection and wrap it with spatial reference
source = osgeo.osr.SpatialReference()
source.SetProjection (ds.GetProjection())

# build the spatial reference for WGS84
target = osgeo.osr.SpatialReference()
target.SetWellKnownGeogCS('WGS84')

# make the transform to go between the CS
transform = osgeo.osr.CoordinateTransformation(source, target)

# take the first point of the map's NEATLINE and transform
p = osgeo.ogr.CreateGeometryFromWkt ('POINT (34277.755142044210515
-376750.367268434783909')
p.Transform (transform)
TypeError                                 Traceback (most recent call
last)
<ipython-input-75-05269f067010> in <module>()
---> 13     p.Transform (transform)
     14     print (p)
     15     # convert to pixel

/usr/lib/python3/dist-packages/osgeo/ogr.py in Transform(self, *args)
   4619         OGRERR_NONE on success or an error code.
   4620         """
-> 4621         return _ogr.Geometry_Transform(self, *args)
   4622 
   4623     def GetSpatialReference(self, *args):

TypeError: in method 'Geometry_Transform', argument 2 of type
'OSRCoordinateTransformationShadow *'

I cannot find anything to help me out at this point -- I am probably not
searching for the right keywords more than anything. Did I mess up the
transform? Did I mess up defining the source projection? I figure these
are the two most likely problems that I have, but I am all ears for any
other problems.

I am worried that all of the suggestions will be to use warp or
transform at the command line, but for complicated reasons using those
tools is undesirable. It is more important that I learn how to do such a
seemingly obvious task in the Python environment.

Thanks in advance for any and all help.

Here are my source and target spatial references:

PROJCS["unnamed",

PROJECTION["PROJCS["unnamed",GEOGCS["NAD83",DATUM["North_American_Datum_1983",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6269"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9108"]],AUTHORITY["EPSG","4269"]],PROJECTION["Albers_Conic_Equal_Area"],PARAMETER["standard_parallel_1",34],PARAMETER["standard_parallel_2",40.5],PARAMETER["latitude_of_center",0],PARAMETER["longitude_of_center",-120],PARAMETER["false_easting",0],PARAMETER["false_northing",-4000000],UNIT["Meter",1]]"]]

GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        TOWGS84[0,0,0,0,0,0,0],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9108"]],
    AUTHORITY["EPSG","4326"]]

-- 
Al Niessner

I have never found the companion that was so companionable as solitude.
- From Walden by Henry David Thoreau

The universe is indifferent, and life is brutal; however, it is man's
choice of behavior that makes them malevolent rather than benevolent.

Some will fall in love with life and drink it from a fountain
That is pouring like an avalanche coming down the mountain.
- From the song Pepper by the Butthole Surfers