Spatial Users Forum

Re: Problem with st_transform

Posted By: ROBERT ULEMAN
Date: Friday, 24 April 2009, at 7:58 p.m.

In Response To: Problem with st_transform (JORGE INFANTE)

Jorge,

Thank you for this question. I'd love to see more traffic on this forum!

Fortunately, this is an easy one. It really helps that you present a complete test case; thank you for that.

You did everything right, but a tiny typographical error tripped you up. In the srtext value you insert into the spatial_ref_sys table, you substitute the WGS 1984 datum for the GRS 1980 datum that actually underlies the POSGAR projections. I don't know the source of your datum specification (perhaps copied from the srtext for srid=4?), but the last two digits in the second spheroid parameter are transposed (_highlights_ inserted): 298.2572235_36_ should be 298.2572235_63_ . When I made that correction, the ST_Transform function worked for me:

> select st_transform(the_geom, 4) from mytable;

(expression) 4 MULTIPOLYGON (((-60.6447115552 -32.95151404, -60.6447043468 -32.9518464663, -60.6443271611 -32.9519433611, -60.644101185 -32.9516707816, -60.6443388958 -32.9514054106, -60.6447115552 -32.95151404)))

1 row(s) retrieved.

As you probably know, the ST_Transform() function can only transform data between spatial reference systems that are defined on the same datum. This is what "incompatible" in the error message refers to: in this case the datums were not the same because their definitions differed in those two transposed digits. Unfortunately, the Spatial DataBlade is not particularly smart about detecting datum differences: it performs a strict character-based comparison, not a numeric one based on knowledge of the meaning of each parameter. Therefore, the character strings must be identical (though I don't know if the comparison is case-sensitive) and EVEN THE PRESENCE OR ABSENCE OF TRAILING ZEROS can make the comparison fail. In this case, the numeric values were actually different, so that needed to be fixed anyway.

A few more comments; these may be arcane to most readers, but I hope some will find them useful.

As you (Jorge) work for a cadastral agency, you are probably well aware that there is a difference between the GRS80 and WGS84 geodetic references: the inverse-flattening spheroid parameters are slightly different. To handle this properly, a real datum transformation is needed; this requires specialized software, and the Informix Spatial DataBlade does not provide this. For many GIS applications, working with data of limited accuracy, the difference is negligible; but for surveying purposes, it is not. If you know what you're doing (and it certainly looks like you do) then you are obviously able to make your own choice whether to accept the discrepancies that may result.

It is certainly not a requirement, but I have found it helpful to choose the SRID as well as the coordinate offsets and units myself, rather than rely on the function SE_CreateSRID(). These days, our internal coordinates are expressed in big (64-bit) integers, so there is almost never a need to fine-tune the offsets and units for individual project areas; a single choice will suffice for all uses of a particular coordinate system. In other words, you can cover the entire earth and still have sub-mm precision. This means that you can tie each SRS you define one-for-one to the associated EPSG coordinate reference, which in turn means that you can use the EPSG SRID (22185, in this case) for your local database SRID as well. That makes it very easy to keep track of what each SRID stands for, and if applied consistently gives SRIDs a fixed meaning across multiple databases.

As to the choice of (X,Y) offsets and units, you chose a range of about 5500 km in X and only 45 km in Y. I'm not sure this is right: Transverse Mercator projections are usually used for regions with a greater north-south extent. In any case, SE_CreateSRID() produced an offset of -537 km and 5,791 km in X and Y, respectively, and an xyunits value of over a billion, for a resolution of less than 1 nm. Perhaps a more widely usable choice would be to follow the published limits (from spatialreference.org) with very wide margins:
Projected Bounds: 5346660.9906, 5660186.1166, 5653339.0094, 7412332.1438
X offset: -1,000,000 (in case you really want to go west to near-zero X-values)
Y offset: -1,000,000 (it doesn't hurt to start so far away from any reasonable coordinates you'll ever encounter)
X/Y units: 1,000,000 (for a resolution of 1 micron)
With this, you still have enough digits to represent positive coordinates in the billions; you have round numbers in the offsets and units, which makes it easy to do back-of-the-envelope calculations and estimates in case you suspect something is wrong; and no matter what kind of buffer operations you do, the internal coordinate representation will never be a limitation. While it may be useful at times to have the "coordinate out of range" error alert you to a bug in a procedure or client, it is generally not a good idea to choose your limits so tight that this becomes a reliable error trap, because if you ever wanted to exceed the official bounds for some algorithmic or computational reason, you would be stuck, with no escape.
With all this, your definition of the projected SRS would use this statement:

INSERT INTO spatial_references (
srid, description, auth_name, auth_srid,
falsex, falsey, xyunits, falsez, zunits, falsem, munits,
srtext
) VALUES (
22185, 'POSGAR 94 / Argentina 5', 'EPSG', 22185,
-1000000, -1000000, 1000000, -1000, 1000, -1000, 1000,
'PROJCS["Gauss_Kruger_Faja5", GEOGCS["GCS_WGS_1984", DATUM["D_WGS_1984",SPHEROID["WGS_1984",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["Degree",0.0174532925199433]],
PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",5500000.0],PARAMETER["False_Northing",0.0],
PARAMETER["Central_Meridian",-60.0],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_of_Origin",-90.0],UNIT["Meter",1.0]]'
);

I hope that I've answered your original question and not made things too confusing with all the other stuff. Don't hesitate to post further comments or questions. If you want to take this off-line, you can contact me at the email address below.

-Robert Uleman
Worldwide Technical Sales
Spatiotemporal Data Management
uleman@us.ibm.com

Messages In This Thread

• Problem with st_transform
  JORGE INFANTE -- Wednesday, 22 April 2009, at 10:54 p.m.
  • Re: Problem with st_transform
    ROBERT ULEMAN -- Friday, 24 April 2009, at 7:58 p.m.

Spatial Users Forum is maintained by Administrator with WebBBS 5.12.