~ubuntu-branches/ubuntu/saucy/merkaartor/saucy

« back to all changes in this revision

Viewing changes to include/ggl/projections/proj/gnom.hpp

  • Committer: Bazaar Package Importer
  • Author(s): Bernd Zeimetz
  • Date: 2009-09-13 00:52:12 UTC
  • mto: (1.2.7 upstream) (0.1.3 upstream) (3.1.7 sid)
  • mto: This revision was merged to the branch mainline in revision 10.
  • Revision ID: james.westby@ubuntu.com-20090913005212-pjecal8zxm07x0fj
ImportĀ upstreamĀ versionĀ 0.14+svnfixes~20090912

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#ifndef GGL_PROJECTIONS_GNOM_HPP
 
2
#define GGL_PROJECTIONS_GNOM_HPP
 
3
 
 
4
// Generic Geometry Library - projections (based on PROJ4)
 
5
// This file is automatically generated. DO NOT EDIT.
 
6
 
 
7
// Copyright Barend Gehrels (1995-2009), Geodan Holding B.V. Amsterdam, the Netherlands.
 
8
// Copyright Bruno Lalande (2008-2009)
 
9
// Use, modification and distribution is subject to the Boost Software License,
 
10
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
 
11
// http://www.boost.org/LICENSE_1_0.txt)
 
12
 
 
13
// This file is converted from PROJ4, http://trac.osgeo.org/proj
 
14
// PROJ4 is originally written by Gerald Evenden (then of the USGS)
 
15
// PROJ4 is maintained by Frank Warmerdam
 
16
// PROJ4 is converted to Geometry Library by Barend Gehrels (Geodan, Amsterdam)
 
17
 
 
18
// Original copyright notice:
 
19
 
 
20
// Permission is hereby granted, free of charge, to any person obtaining a
 
21
// copy of this software and associated documentation files (the "Software"),
 
22
// to deal in the Software without restriction, including without limitation
 
23
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
 
24
// and/or sell copies of the Software, and to permit persons to whom the
 
25
// Software is furnished to do so, subject to the following conditions:
 
26
 
 
27
// The above copyright notice and this permission notice shall be included
 
28
// in all copies or substantial portions of the Software.
 
29
 
 
30
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
 
31
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 
32
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
 
33
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 
34
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
 
35
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
 
36
// DEALINGS IN THE SOFTWARE.
 
37
 
 
38
#include <boost/math/special_functions/hypot.hpp>
 
39
 
 
40
#include <ggl/projections/impl/base_static.hpp>
 
41
#include <ggl/projections/impl/base_dynamic.hpp>
 
42
#include <ggl/projections/impl/projects.hpp>
 
43
#include <ggl/projections/impl/factory_entry.hpp>
 
44
 
 
45
namespace ggl { namespace projection
 
46
{
 
47
    #ifndef DOXYGEN_NO_DETAIL
 
48
    namespace detail { namespace gnom{ 
 
49
            static const double EPS10 = 1.e-10;
 
50
            static const int N_POLE = 0;
 
51
            static const int S_POLE = 1;
 
52
            static const int EQUIT = 2;
 
53
            static const int OBLIQ = 3;
 
54
 
 
55
            struct par_gnom
 
56
            {
 
57
                double sinph0;
 
58
                double cosph0;
 
59
                int  mode;
 
60
            };
 
61
 
 
62
            // template class, using CRTP to implement forward/inverse
 
63
            template <typename Geographic, typename Cartesian, typename Parameters>
 
64
            struct base_gnom_spheroid : public base_t_fi<base_gnom_spheroid<Geographic, Cartesian, Parameters>,
 
65
                     Geographic, Cartesian, Parameters>
 
66
            {
 
67
 
 
68
                 typedef double geographic_type;
 
69
                 typedef double cartesian_type;
 
70
 
 
71
                par_gnom m_proj_parm;
 
72
 
 
73
                inline base_gnom_spheroid(const Parameters& par)
 
74
                    : base_t_fi<base_gnom_spheroid<Geographic, Cartesian, Parameters>,
 
75
                     Geographic, Cartesian, Parameters>(*this, par) {}
 
76
 
 
77
                inline void fwd(geographic_type& lp_lon, geographic_type& lp_lat, cartesian_type& xy_x, cartesian_type& xy_y) const
 
78
                {
 
79
                        double  coslam, cosphi, sinphi;
 
80
                
 
81
                        sinphi = sin(lp_lat);
 
82
                        cosphi = cos(lp_lat);
 
83
                        coslam = cos(lp_lon);
 
84
                        switch (this->m_proj_parm.mode) {
 
85
                        case EQUIT:
 
86
                                xy_y = cosphi * coslam;
 
87
                                break;
 
88
                        case OBLIQ:
 
89
                                xy_y = this->m_proj_parm.sinph0 * sinphi + this->m_proj_parm.cosph0 * cosphi * coslam;
 
90
                                break;
 
91
                        case S_POLE:
 
92
                                xy_y = - sinphi;
 
93
                                break;
 
94
                        case N_POLE:
 
95
                                xy_y = sinphi;
 
96
                                break;
 
97
                        }
 
98
                        if (xy_y <= EPS10) throw proj_exception();;
 
99
                        xy_x = (xy_y = 1. / xy_y) * cosphi * sin(lp_lon);
 
100
                        switch (this->m_proj_parm.mode) {
 
101
                        case EQUIT:
 
102
                                xy_y *= sinphi;
 
103
                                break;
 
104
                        case OBLIQ:
 
105
                                xy_y *= this->m_proj_parm.cosph0 * sinphi - this->m_proj_parm.sinph0 * cosphi * coslam;
 
106
                                break;
 
107
                        case N_POLE:
 
108
                                coslam = - coslam;
 
109
                        case S_POLE:
 
110
                                xy_y *= cosphi * coslam;
 
111
                                break;
 
112
                        }
 
113
                }
 
114
 
 
115
                inline void inv(cartesian_type& xy_x, cartesian_type& xy_y, geographic_type& lp_lon, geographic_type& lp_lat) const
 
116
                {
 
117
                        double  rh, cosz, sinz;
 
118
                
 
119
                        rh = boost::math::hypot(xy_x, xy_y);
 
120
                        sinz = sin(lp_lat = atan(rh));
 
121
                        cosz = sqrt(1. - sinz * sinz);
 
122
                        if (fabs(rh) <= EPS10) {
 
123
                                lp_lat = this->m_par.phi0;
 
124
                                lp_lon = 0.;
 
125
                        } else {
 
126
                                switch (this->m_proj_parm.mode) {
 
127
                                case OBLIQ:
 
128
                                        lp_lat = cosz * this->m_proj_parm.sinph0 + xy_y * sinz * this->m_proj_parm.cosph0 / rh;
 
129
                                        if (fabs(lp_lat) >= 1.)
 
130
                                                lp_lat = lp_lat > 0. ? HALFPI : - HALFPI;
 
131
                                        else
 
132
                                                lp_lat = asin(lp_lat);
 
133
                                        xy_y = (cosz - this->m_proj_parm.sinph0 * sin(lp_lat)) * rh;
 
134
                                        xy_x *= sinz * this->m_proj_parm.cosph0;
 
135
                                        break;
 
136
                                case EQUIT:
 
137
                                        lp_lat = xy_y * sinz / rh;
 
138
                                        if (fabs(lp_lat) >= 1.)
 
139
                                                lp_lat = lp_lat > 0. ? HALFPI : - HALFPI;
 
140
                                        else
 
141
                                                lp_lat = asin(lp_lat);
 
142
                                        xy_y = cosz * rh;
 
143
                                        xy_x *= sinz;
 
144
                                        break;
 
145
                                case S_POLE:
 
146
                                        lp_lat -= HALFPI;
 
147
                                        break;
 
148
                                case N_POLE:
 
149
                                        lp_lat = HALFPI - lp_lat;
 
150
                                        xy_y = -xy_y;
 
151
                                        break;
 
152
                                }
 
153
                                lp_lon = atan2(xy_x, xy_y);
 
154
                        }
 
155
                }
 
156
            };
 
157
 
 
158
            // Gnomonic
 
159
            template <typename Parameters>
 
160
            void setup_gnom(Parameters& par, par_gnom& proj_parm)
 
161
            {
 
162
                if (fabs(fabs(par.phi0) - HALFPI) < EPS10)
 
163
                        proj_parm.mode = par.phi0 < 0. ? S_POLE : N_POLE;
 
164
                else if (fabs(par.phi0) < EPS10)
 
165
                        proj_parm.mode = EQUIT;
 
166
                else {
 
167
                        proj_parm.mode = OBLIQ;
 
168
                        proj_parm.sinph0 = sin(par.phi0);
 
169
                        proj_parm.cosph0 = cos(par.phi0);
 
170
                }
 
171
                // par.inv = s_inverse;
 
172
                // par.fwd = s_forward;
 
173
                par.es = 0.;
 
174
            }
 
175
 
 
176
        }} // namespace detail::gnom
 
177
    #endif // doxygen 
 
178
 
 
179
    /*!
 
180
        \brief Gnomonic projection
 
181
        \ingroup projections
 
182
        \tparam Geographic latlong point type
 
183
        \tparam Cartesian xy point type
 
184
        \tparam Parameters parameter type
 
185
        \par Projection characteristics
 
186
         - Azimuthal
 
187
         - Spheroid
 
188
        \par Example
 
189
        \image html ex_gnom.gif
 
190
    */
 
191
    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
 
192
    struct gnom_spheroid : public detail::gnom::base_gnom_spheroid<Geographic, Cartesian, Parameters>
 
193
    {
 
194
        inline gnom_spheroid(const Parameters& par) : detail::gnom::base_gnom_spheroid<Geographic, Cartesian, Parameters>(par)
 
195
        {
 
196
            detail::gnom::setup_gnom(this->m_par, this->m_proj_parm);
 
197
        }
 
198
    };
 
199
 
 
200
    #ifndef DOXYGEN_NO_DETAIL
 
201
    namespace detail
 
202
    {
 
203
 
 
204
        // Factory entry(s)
 
205
        template <typename Geographic, typename Cartesian, typename Parameters>
 
206
        class gnom_entry : public detail::factory_entry<Geographic, Cartesian, Parameters>
 
207
        {
 
208
            public :
 
209
                virtual projection<Geographic, Cartesian>* create_new(const Parameters& par) const
 
210
                {
 
211
                    return new base_v_fi<gnom_spheroid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
 
212
                }
 
213
        };
 
214
 
 
215
        template <typename Geographic, typename Cartesian, typename Parameters>
 
216
        inline void gnom_init(detail::base_factory<Geographic, Cartesian, Parameters>& factory)
 
217
        {
 
218
            factory.add_to_factory("gnom", new gnom_entry<Geographic, Cartesian, Parameters>);
 
219
        }
 
220
 
 
221
    } // namespace detail 
 
222
    #endif // doxygen
 
223
 
 
224
}} // namespace ggl::projection
 
225
 
 
226
#endif // GGL_PROJECTIONS_GNOM_HPP
 
227