# Dozier's TM method---my summary

11 messages
Open this post in threaded view
|

## Dozier's TM method---my summary

 I think the basic idea presented by Dozier sounds feasible but his execution falls short of the goal. Although the the method of computing elliptic integral of the second kind with complex argument may be OK it has troubles with large arguments when evaluating the Jacobian Zeta function.  The second, and I believe the more serious problem, is the Newton-Raphson method employed.  He has expanded the basic real function and applied it to a complex variable.  I am not sure that this is appropriate and searching through the net and all leads me to the conclusion that we are getting into deep water when dealing with the complex plane.  Also, it looks like we are also dealing with multiple roots, especially when longitude exceeds a certain value (suggested to be (pi/2)*(1-k))---a factor not addressed in Dozier's solution. Note that the cusp (poorly displayed in the previously mention gif url) appears to be the beginning of the multiple root solution.  I say poorly displayed as the gradation of the equitorial parallel should *smoothly* begin a swing to the north OR (importantly) to the south.  Selecting the north or south root becomes a practical problem for a projection program and is a problem with any method dealing with the comprehensive TM projection. Because I have no training and no experience in working with complex variable problems and have failed to find any practical material related to alternate methods to compute elliptic integrals with complex arguments and, more importantly, a Newton-Raphson routine for determining roots of complex functions, I have decide to suspend again any activity on the Dozier method. One has to know when to throw in the towel.  ;-) -- Jerry and the low-riders: Daisy Mae and Joshua "Cogito cogito ergo cogito sum"    Ambrose Bierce, The Devil's Dictionary _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Dozier's TM method---my summary

Open this post in threaded view
|

## Re: Re: Dozier's TM method---my summary

 On Thursday 29 June 2006 3:30 pm, [hidden email] wrote:         ... > > Note that the cusp (poorly displayed in the previously mention gif url) > > appears to be the beginning of the multiple root solution.  I say poorly > > displayed as the gradation of the equitorial parallel should *smoothly* > > begin > > a swing to the north OR (importantly) to the south.  > > The cusp is poorly observed, not poorly displayed. The illustration is > exact to sub-pixel resolution and it matches Lee's illustration precisely > despite using a different method of calculation. It seems you do not know a > continuously differentiable curve when you see one. My goodness, such a cranky old man!   I know you don't like me but I have tried to be civil.  Can we bury the hatchet sometime.  I am getting too old for this bullshit. -- Jerry and the low-riders: Daisy Mae and Joshua "Cogito cogito ergo cogito sum"    Ambrose Bierce, The Devil's Dictionary _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Re: Dozier's TM method---my summary

Open this post in threaded view
|

## Re: Re: Dozier's TM method---my summary

 > From: "Gerald I. Evenden" > Stating it another way, my surrender was to the Dozier article that while the > method sounds interesting the execution falls short.  The code had several > errors and some methodology is in question.  I have spent an additional week > rooting out a complex Bulirsch routine from the old Bell lab site and spent > several days relating Jacobi's form to Legendre in order to test against > Abramowitz tables and GSL software for at least real arguments---they now > agree.  But Dozier's code does not agree to better than the third or forth > place (I am not saying Dozier's code is wrong yet). > > And I still have no handle on a replacement for Newton-Raphson. What surprises me is that I gave an alternative twice, once in a posting on this list dated 2006-06-13, once in a private email dated 2006-06-23, including a web reference. It's about ACM TOMS algorithm 365 from H. Bach. Fortran code can be found at netlib.org and one can buy a copy of the peer reviewed TOMS article at ACM (acm.org). I am not happy with the method because it is horribly slow, but it is a very robust method, at least if one knows how to tweak its control parameters. As for a replacement of complex Jacobi elliptic functions and a complex elliptic integral: it can be done without too much effort. I provided in the same manner (list posting and private email) references. I do admit that working with complex elliptics is no trivial endeavor without knowledge or experience. _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Re: Dozier's TM method---my summary

 On Wednesday 05 July 2006 8:31 am, Oscar van Vlijmen wrote: > > From: "Gerald I. Evenden" > > > > Stating it another way, my surrender was to the Dozier article that while > > the method sounds interesting the execution falls short.  The code had > > several errors and some methodology is in question.  I have spent an > > additional week rooting out a complex Bulirsch routine from the old Bell > > lab site and spent several days relating Jacobi's form to Legendre in > > order to test against Abramowitz tables and GSL software for at least > > real arguments---they now agree.  But Dozier's code does not agree to > > better than the third or forth place (I am not saying Dozier's code is > > wrong yet). > > > > And I still have no handle on a replacement for Newton-Raphson. > > What surprises me is that I gave an alternative twice, once in a posting on > this list dated 2006-06-13, once in a private email dated 2006-06-23, > including a web reference. > It's about ACM TOMS algorithm 365 from H. Bach. Fortran code can be found > at netlib.org and one can buy a copy of the peer reviewed TOMS article at > ACM (acm.org). > I am not happy with the method because it is horribly slow, but it is a > very robust method, at least if one knows how to tweak its control > parameters. > > As for a replacement of complex Jacobi elliptic functions and a complex > elliptic integral: it can be done without too much effort. I provided in > the same manner (list posting and private email) references. > I do admit that working with complex elliptics is no trivial endeavor > without knowledge or experience. The Bulirsch code seems to work and compares in the real case (x+0i) with the GSL library *and* Abramowich tables--- after a lot of footwork dealing with figuring out how to match up the (&*^(&^ parameters to each.  The input parameters were what was really blowing my mind. I ran into another, perhaps the same source, FORTRAN source but found out I did not have a g77 compiler (I had recently  changed versions of Linux) and solving that problem introduced more headaches than I wanted at the moment. In my little bench program where I compared GSL, Bulirsch and Dozier I ran into consistant disagreement with Dozier code in the 3rd or 4th decimal place.  At the moment, this does not make sense because Dozier's code is at a significant point in Dozier's second phase of the forward solution I do not understand how it could be incorrect and yet yield correct answers in many areas for the overall projection. The following is a short example of a test run: e^2 used: 0.006722670022 x + y*i input: 0.1 0 Dozier e2ci: 0.09999998496533813  0*i Bulirsch elco2 ret: OK elco2: 0.0999988817825081  0*i gsl_sf_ellint_E_e x: 0.09999888178250819 --- err: 4.44104e-17 x + y*i input: 0.2 0 Dozier e2ci: 0.19999988044241  0*i Bulirsch elco2 ret: OK elco2: 0.1999911075211161  0*i gsl_sf_ellint_E_e x: 0.1999911075211167 --- err: 8.88297e-17 x + y*i input: 0.5 0 Dozier e2ci: 0.4999982088519777  0*i Bulirsch elco2 ret: OK elco2: 0.4998667513591754  0*i gsl_sf_ellint_E_e x: 0.4998667513591855 --- err: 2.22222e-16 x + y*i input: 0.3 0.3 Dozier e2ci: 0.3000007837129038  0.2999991577227006*i Bulirsch elco2 ret: OK elco2: 0.3000583008002545  0.2999373655430532*i In the last case, GSL not run because y != 0. That gives you an idea where I now stand with the complex elliptical integral of the second kind (incomplete). The nasty part is that even if this is solved, the problem with the first phase and the root finding looms as an even bigger issue. Sigh.  :-( > _______________________________________________ > Proj mailing list > [hidden email] > http://lists.maptools.org/mailman/listinfo/proj-- Jerry and the low-riders: Daisy Mae and Joshua "Cogito cogito ergo cogito sum"    Ambrose Bierce, The Devil's Dictionary _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Re: Dozier's TM method---my summary

 In reply to this post by OvV_HN On Wednesday 05 July 2006 8:31 am, Oscar van Vlijmen wrote:         ... > What surprises me is that I gave an alternative twice, once in a posting on > this list dated 2006-06-13, once in a private email dated 2006-06-23, > including a web reference. > It's about ACM TOMS algorithm 365 from H. Bach. Fortran code can be found > at netlib.org and one can buy a copy of the peer reviewed TOMS article at > ACM (acm.org). > I am not happy with the method because it is horribly slow, but it is a > very robust method, at least if one knows how to tweak its control > parameters. Not sure what you meant by slow but the example clocked in at 2ms real time and did not show on the millisecond user/system clocks. Sorta looks like a steepest descent routine with no analysis from derivatives so it is likely to be sluggish. So much for Dozier's idea of a speedy method.  ;-) I stupidly missed all the control in the file so had to do about three whacks at getting it to compile.  Once the junk was removed it went fine. -- Jerry and the low-riders: Daisy Mae and Joshua "Cogito cogito ergo cogito sum"    Ambrose Bierce, The Devil's Dictionary _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Using libproj with array data

 Hi, I hope this is not a FAQ or covered in the existing docs; I've   checked and apologize if I missed something obvious. I have a data server which is used to access a number of files of   georeferenced data stored in two dimentional arrays. For each array,   I know the longitude and latitude values at the corners as well as   the projection and datum. Is there a way I can use libproj to map   long/lat locations to indices within the arrays? Thanks, James -- James Gallagher                jgallagher at opendap.org OPeNDAP, Inc                   406.723.8663 _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj
Open this post in threaded view
|

## Re: Using libproj with array data

 James, What you're trying to do sounds like re-projecting a grid. (A grid just being a two dimensional array of values). I have written code to do this previously. To project a grid you basically have to perform a forward projection of some points around the edge of the input grid and use these to compute the bounding rectangle in the new coordinate system. (Remember some projections will distort the shape of the grid). You then create an output grid in the new coordinate system with the appropriate dimensions. Finally you have to go through each cell in the output grid and use reverse projection to find the corresponding cell in the unprojected grid. Use grid interpolation to pick up the value you need to put in the output grid cell. There's probably already code in the public domain which uses PROJ to do this. Anthony. --- James Gallagher <[hidden email]> wrote: > Hi, I hope this is not a FAQ or covered in the existing docs; I've   > checked and apologize if I missed something obvious. > > I have a data server which is used to access a number of files of   > georeferenced data stored in two dimentional arrays. For each array, > > I know the longitude and latitude values at the corners as well as   > the projection and datum. Is there a way I can use libproj to map   > long/lat locations to indices within the arrays? > > Thanks, > James > -- > James Gallagher                jgallagher at opendap.org > OPeNDAP, Inc                   406.723.8663 > > _______________________________________________ > Proj mailing list > [hidden email] > http://lists.maptools.org/mailman/listinfo/proj> __________________________________________________ Do You Yahoo!? Tired of spam?  Yahoo! Mail has the best spam protection around http://mail.yahoo.com  _______________________________________________ Proj mailing list [hidden email] http://lists.maptools.org/mailman/listinfo/proj