ocea - Incorrect answers

[proj4-bot]

Reported by fwillett on 10 May 2012 17:41 UTC
System:
OS - Windows 7 Ultimate, SP1, 64 bit
CPU - Intel Core2 Quad CPU @ 2.83 GHz
Memory - 1.86 GB
Running in Oracle VirtualBox 4.16 r74713

Installation:
proj version - Rel. 4.8.0, 6 March 2012
Microsoft Visual Studio 2010 Professional Version 10.0.40219.1 SP1Rel
Build commands:
"C:\Program Files (x86)\Microsoft Visual Studio 10.0\VC\bin\vcvars32.bat"
nmake /f makefile.vc
nmake /f makefile.vc install-all

Note: This builds a 32 bit version of the programs and libraries.

In the following the OS prompt is "$" and the proj prompt is "%"

ocea - Incorrect answers
Using original code:
$proj +proj=ocea +ellps=sphere +lat_1=-45 +lat_2=45 +lon_1=-90 +lon_2=90
%0 0
0.00 0.00
%0.0001 0
50093014.90 -50093014.90
%0 0.0001
50093014.90 50093014.90

Appied the following patch to PJ_ocea.c
// http://osgeo-org.1560.n6.nabble.com/Oblique-Cyl-Eq-Area-td3851618.html
// Gerald I. Evenden - Jul 07, 2005
// Alas, the solution to the problem if you wish to fix your own copy of the source is to
// remove the reference to the earth's radius in PJ_ocea.c. The first two executable
// statements just after the entry ENTRY0(ocea) read
//
// P->rok = P->a / P->k0;
// P->rtk = P->a * P->k0;
//
// Change the "P->a" to "1." on both lines. Simplify the second line if you like. I suspect
// that this error has persisted over a long line of releases and versions and points to the
// folly of using a unity earth radius as a test value (see p. 279, PP 1395). Also shows
// how often ocea is used. This is such a dumb error. Only one projection references the
// earth's radius (or major axis).
//
// P->rok = P->a / P->k0;
// P->rtk = P->a * P->k0;
P->rok = 1. / P->k0;
P->rtk = 1. * P->k0;

New results which appear to be correct:
$proj +proj=ocea +ellps=sphere +lat_1=-45 +lat_2=45 +lon_1=-90 +lon_2=90
%0 0
0.00 0.00
%0.0001 0.0001
15.73 0.00
%0.0001 0
7.86 -7.86
%0 0.0001
7.86 7.86

This modification also seems to fix the alternate way to specify the projection with
+lat_0, +lonc, +alpha, +x_0, +y_0.

Migrated-From: https://trac.osgeo.org/proj/ticket/166

[kbevers]

I wanted to apply the fix mentioned above, but could not reproduce the initial bug. Instead I get the following output when projecting the origin:

$ ./src/proj +proj=ocea +ellps=sphere +lat_1=-45 +lat_2=45 +lon_1=-90 +lon_2=90
0 0
127515997886954.58  0.00

It wasn't immediately clear to me when this new error was introduced. With git bisect I've managed to narrow it down to using DEG_TO_RAD instead of a using the literal number in dmstor.c. This change was introduced in 86b2155 by @cffk. It's far too close to my bedtime for me to understand why that should cause a problem in PJ_ocea.. That's a problem for another day...