[GRASSLIST:1481] Re: Bug in s.proj?

[David Orme]

Hi,

OK - I've followed Hamish's suggestion and used cs2cs from PROJ.4 (I'll 
get to Grass 5.3 soon) to check the location of a test point (0 m E, 
7332956.368 m N) from a Behrmann projection back into geographic 
coordinates:

[doibook:/usr/local/bin] dorme% ./cs2cs +proj=cea +lat_ts=30 
+datum=WGS84 +to +proj=latlong +datum=WGS84
0 7332956.368
0dE     86d18'12.382"N 0.000

This gives me the same geographic coordinate as I got from my original 
use of s.proj in Grass 5.0.2. To test ArcMAP again, I created a 
shapefile with a Behrmann projection containing a single point at the 
test location and used the project wizard in ArcToolbox to change to a 
Geographic projection (wgs84.prj). The resulting point has the 
following geographic coordinates:

0dE 87d 7' 58.42" N

The output of g.projinfo for my Behrmann location is below. I'm using 
Mac OS X (version 10.2.8) and the following version of cs2cs:

[doibook:/usr/local/bin] dorme% ./cs2cs
Rel. 4.4.7, 31 March 2003
usage: cs2cs [ -eEfIlrstvwW [args] ] [ +opts[=arg] ]
                    [+to [+opts[=arg] [ files ]

David.

-----------------------------------------------------------
PROJ_INFO file:
name:    Equal Area Cylindrical
datum:   wgs84
dx:      0.000000
dy:      0.000000
dz:      0.000000
proj:    cea
ellps:   wgs84
a:       6378137.0000000000
es:      0.0066943800
f:       298.2572235630
lat_0:   0.0000000000
lat_ts:  30.0000000000
lon_0:   0.0000000000
-----------------------------------------------------------
PROJ_UNITS file:
unit:    meter
units:   meters
meters:  1.000000
-----------------------------------------------------------

On Friday, October 10, 2003, at 10:50  am, Hamish wrote:

>> I e-mailed yesterday about a mismatch between Arcmap and Grass 5.0.2
>> transformations of sites from a global cylindrical equal area location
>> to a geographic. I've been playing around some more and it looks to me
>> as if either:
>>
>> 1) s.proj has a bug in it,
>> 2) arcmap has a bug in the transformation of sites or,
>> 3) I've missed something elementary.
>>
>> I created a point shapefile in arcmap of the centre points of the grid
>> cells for the CEA location and read this into grass using s.in.shape.
>> I've checked that the values in the sites file are what I was
>> expecting. I used s.proj to project these sites into a geographic
>> (wgs84) location. I also used arcmap to project the original behrmann
>> points shapefile into wgs84 and then imported this into my grass wgs84
>> location using s.in.shape. The two sets of sites are not in the same
>> place - in extreme latitudes the discrepancy is about half a degree. I
>> won't attach the files to preserve bandwidth but if anyone wants them,
>> let me know.
>>
>> Any suggestions?
>
>
> 5.0.x does not do datum transforms.
>
> Did you try GRASS 5.3? Is it still broken there?
> Half a degree is pretty huge though.
> Try with 5.3's m.proj2?
>   echo "easting northing" | m.proj2 -o
>
> Try with PROJ.4's cs2cs? (what grass 5.3 uses more or less)
>   http://www.remotesensing.org/proj
>
>
>
> If it persists after trying 5.3 can you file a bug at:
> http://grass.ibiblio.org/bugtracking/bugreport.html
> ?
> please include some specific test points with the results from both
> GISs, output of g.projinfo, and any info you have on the projection
> parameters.
>
> If you can't try 5.3, just post a bug report with the test point and
> someone else can check.
>
>
>
> thanks,
> Hamish
>