# Problems with Pittman geodesic??

8 messages
Open this post in threaded view
|

## Problems with Pittman geodesic??

 While checking out the accuracy of Vincenty vs. Pittman I was fairly consistently getting agreement to micron or better. BUT not always!!! I would appreciate anyone who has a Pittman geodesic routine to test the following inverse (WGS84 ellipsoid): point 1 at 0 lat, 0 lon and point 2 at 45N 90E. With my version of Pittman I get a distance of 9993541.5348708 AND I get the same distance while changing the longitude of the second point from about 89.6 to a hair over 90. At 0,0/45,89.5: Pittman: 9970963.0100082 Vincenty:9970963.01000801 At 0,0/45,89.6 Pittman: 9993541.5348708 Vincenty:9978847.65947167 Pittman remained constant in this longitude interval At 0,0/45,90 Pittman:  9993541.5348708 Vincenty:10010386.3610382 Between 90.0000001 and 90.000001 Pittman finally came back into near agreement with Vincenty. Note: I did nothing to the Pittman FORTRAN loaned from Mugnier other than compile it with gfortran and link to a C driver. -- The whole religious complexion of the modern world is due to the absence from Jerusalem of a lunatic asylum. -- Havelock Ellis (1859-1939) British psychologist _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Problems with Pittman geodesic??

 Gerald wrote: > While checking out the accuracy of Vincenty vs. Pittman > I was fairly consistently getting agreement to micron or > better. > > BUT not always!!! > > I would appreciate anyone who has a Pittman geodesic routine > to test the following inverse (WGS84 ellipsoid): > > point 1 at 0 lat, 0 lon and point 2 at 45N 90E. > > With my version of Pittman I get a distance of 9993541.5348708 > AND I get the same distance while changing the longitude of the > second point from about 89.6 to a hair over 90. > > At 0,0/45,89.5: > Pittman: 9970963.0100082 > Vincenty:9970963.01000801 > > At 0,0/45,89.6 > Pittman: 9993541.5348708 > Vincenty:9978847.65947167 > > Pittman remained constant in this longitude interval > > At 0,0/45,90 > Pittman:  9993541.5348708 > Vincenty:10010386.3610382 > > Between 90.0000001 and 90.000001 Pittman finally came back > into near agreement with Vincenty. > > Note: I did nothing to the Pittman FORTRAN loaned from Mugnier > other than compile it with gfortran and link to a C driver. I would guess that some internal variable is defined as a REAL (float) instead of a DOUBLE PRECISION (double), and so the effect you are seeing is quantization due to variable precision. Hamish       _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Problems with Pittman geodesic??

 On Sunday 18 January 2009 11:00:36 pm Hamish wrote: > Gerald wrote: > > While checking out the accuracy of Vincenty vs. Pittman > > I was fairly consistently getting agreement to micron or > > better. > > > > BUT not always!!! > > > > I would appreciate anyone who has a Pittman geodesic routine > > to test the following inverse (WGS84 ellipsoid): > > > > point 1 at 0 lat, 0 lon and point 2 at 45N 90E. > > > > With my version of Pittman I get a distance of 9993541.5348708 > > AND I get the same distance while changing the longitude of the > > second point from about 89.6 to a hair over 90. > > > > At 0,0/45,89.5: > > Pittman: 9970963.0100082 > > Vincenty:9970963.01000801 > > > > At 0,0/45,89.6 > > Pittman: 9993541.5348708 > > Vincenty:9978847.65947167 > > > > Pittman remained constant in this longitude interval > > > > At 0,0/45,90 > > Pittman:  9993541.5348708 > > Vincenty:10010386.3610382 > > > > Between 90.0000001 and 90.000001 Pittman finally came back > > into near agreement with Vincenty. > > > > Note: I did nothing to the Pittman FORTRAN loaned from Mugnier > > other than compile it with gfortran and link to a C driver. > > I would guess that some internal variable is defined as a REAL (float) > instead of a DOUBLE PRECISION (double), and so the effect you are seeing > is quantization due to variable precision. Both replies, so far seemed to think I was trying to debug the Pittman procedures. I am not interested in fixing the Pittman procedures.  I am asking for anyone with the Pittman routines to check if their version has the same problem mine has.  My only interest in Pittman is as a comparison to the Vincenty procedures because Pittman has been touted as the standard for geodesic calculations. Thank-you. -- The whole religious complexion of the modern world is due to the absence from Jerusalem of a lunatic asylum. -- Havelock Ellis (1859-1939) British psychologist _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Problems with Pittman geodesic??

 In reply to this post by Gerald I. Evenden Some time ago I referred to 2 articles by Saito. He computes the geodesic problems by numerical integration. There are some errors in the formulae, but his second article gives correct results for 'normal' values. Saito attends much detail to borderline situations, but there are probably some errors in his descriptions. Anyhow, he gives a couple of test points which everybody with geodesics code definitively should try to replicate. You might be surprised! Saito article 2 examples. Using the formulas and the procedure described in the article, with an equatorial radius of 6378388 meters and an inverse flattening of 297, the following results were obtained. The computation was made in double precision arithmetic - through twenty significant digits. Line lat1                lat2               delta lon 1. 37d 19' 54.95367"   26d 07' 42.83946"  041d 28' 35.50729" 2. 35  16  11.24862    67  22  14.77638   137  47  28.31435 3. 01  00  00.00000   -00  59  53.83076   179  17  48.02997 4. 01  00  00.00000    01  01  15.18952   179  46  17.84244 5. 41  41  45.88000    41  41  46.20000   000  00  00.56000 6. 30  00  00.00000    37  53  32.46584   116  19  16.68843 7. 30  19  54.95367   -30  11  50.15681   179  58  17.84244 8. 00  39  49.12586   -00  45  14.13112   179  58  17.84244 9. 00  00  54.95367    00  00  42.83946   179  28  17.84244 Line  dist.           azim.                back azim. 1.  4085966.7026 m   095d 27' 59.630899"  118d 05' 58.961609" 2.  8084823.8383     015  44  23.748498   114  55  39.921473 3. 19959999.9998     088  59  59.998971   091  00  06.118356 4. 19780006.5588     004  59  59.999957   174  59  59.884800 5.       16.2839751  052  40  39.390671   052  40  39.763172 6. 10002499.9999     045  00  00.000004   129  08  12.326009 7. 19989590.5480     002  23  52.108130   177  36  19.670109 8. 19994529.4446     177  39  39.010104   002  20  21.153036 9. 19977290.7711     054  08  27.731619   125  51  32.272327 From: THE COMPUTATION OF LONG GEODESICS ON THE ELLIPSOID THROUGH GAUSSlAN QUADRATURE Tsutomu SAITO Construction College, Kodaira-shi Tokyo, 187 Japan. Oscar van Vlijmen ----- Original Message ----- From: "Gerald I. Evenden" <[hidden email]> To: "PROJ.4 and general Projections Discussions" <[hidden email]> Sent: Monday, January 19, 2009 2:57 AM Subject: [Proj] Problems with Pittman geodesic?? > While checking out the accuracy of Vincenty vs. Pittman > I was fairly consistently getting agreement to micron or > better. > > BUT not always!!! > > I would appreciate anyone who has a Pittman geodesic routine > to test the following inverse (WGS84 ellipsoid): > ..... _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Problems with Pittman geodesic??

 On Monday 19 January 2009 1:16:38 pm OvV_HN wrote: > Some time ago I referred to 2 articles by Saito. > He computes the geodesic problems by numerical integration. There are some > errors in the formulae, but his second article gives correct results for > 'normal' values. Saito attends much detail to borderline situations, but > there are probably some errors in his descriptions. > Anyhow, he gives a couple of test points which everybody with geodesics > code definitively should try to replicate. You might be surprised! I did not find any surprises here because where problems do show up are for the near apodal points which Vincenty (and others) usually flunk. I have more problems with failures with Pittman and some non-apodal points. > Saito article 2 examples. > > Using the formulas and the procedure described in the article, with an > equatorial radius of 6378388 meters and an inverse flattening of 297, the > following results were obtained. The computation was made in double > precision arithmetic - through twenty significant digits. "double precision" --- "through twenty significant digits." ???  What kind of hardware is being used?  "double" on most 64bit systems is about 15 digits +-.  With so much significance, why are the following tables showing 12 at best? > Line lat1                lat2               delta lon > 1. 37d 19' 54.95367"   26d 07' 42.83946"  041d 28' 35.50729" > 2. 35  16  11.24862    67  22  14.77638   137  47  28.31435 > 3. 01  00  00.00000   -00  59  53.83076   179  17  48.02997 > 4. 01  00  00.00000    01  01  15.18952   179  46  17.84244 > 5. 41  41  45.88000    41  41  46.20000   000  00  00.56000 > 6. 30  00  00.00000    37  53  32.46584   116  19  16.68843 > 7. 30  19  54.95367   -30  11  50.15681   179  58  17.84244 > 8. 00  39  49.12586   -00  45  14.13112   179  58  17.84244 > 9. 00  00  54.95367    00  00  42.83946   179  28  17.84244 > > Line  dist.           azim.                back azim. > 1.  4085966.7026 m   095d 27' 59.630899"  118d 05' 58.961609" > 2.  8084823.8383     015  44  23.748498   114  55  39.921473 > 3. 19959999.9998     088  59  59.998971   091  00  06.118356 > 4. 19780006.5588     004  59  59.999957   174  59  59.884800 > 5.       16.2839751  052  40  39.390671   052  40  39.763172 > 6. 10002499.9999     045  00  00.000004   129  08  12.326009 > 7. 19989590.5480     002  23  52.108130   177  36  19.670109 > 8. 19994529.4446     177  39  39.010104   002  20  21.153036 > 9. 19977290.7711     054  08  27.731619   125  51  32.272327 > > From: > THE COMPUTATION OF LONG GEODESICS ON THE ELLIPSOID > THROUGH GAUSSlAN QUADRATURE > Tsutomu SAITO Using Gaussian Quadrature is not a big thing as long as one has the integrals to evaluate---one does not need to use gsl to do that. ;-)  I am surprised that he can get away with just using that method---or does he? Lastly, the Saito source is apparently available though "Springer -- something" that charges \$35(presumably US) for a pdf copy of Vincenty that you can get for free from NOAA. -- The whole religious complexion of the modern world is due to the absence from Jerusalem of a lunatic asylum. -- Havelock Ellis (1859-1939) British psychologist _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Problems with Pittman geodesic??

 From: "Gerald I. Evenden" <[hidden email]> Sent: Thursday, January 22, 2009 9:40 PM Subject: Re: [Proj] Problems with Pittman geodesic?? > On Monday 19 January 2009 1:16:38 pm OvV_HN wrote: >> Some time ago I referred to 2 articles by Saito. .... >> Anyhow, he gives a couple of test points which everybody with geodesics >> code definitively should try to replicate. You might be surprised! > I did not find any surprises here because where problems do show up are > for > the near apodal points which Vincenty (and others) usually flunk. Nice. Indeed my implementation of NOAAs Vincenty went wrong with the antipodal examples 7, 8, 9. As did the earlier mentioned method of Xue-Lian Zhang. >> Saito article 2 examples. ... >> following results were obtained. The computation was made in double >> precision arithmetic - through twenty significant digits. > "double precision" --- "through twenty significant digits." ???  What kind > of > hardware is being used? The article is from 1978; only Fortran is mentioned. > Using Gaussian Quadrature is not a big thing as long as one has the > integrals > to evaluate---one does not need to use gsl to do that. ;-) Then you probably don't need GNU GSL for the Mayr-Tobler projection! Oscar van Vlijmen _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj