Dogcows Code - chaz/yoink/blob - src/cml/mathlib/coord_conversion.h

   1 /* -*- C++ -*- ------------------------------------------------------------
   2
   3 Copyright (c) 2007 Jesse Anders and Demian Nave http://cmldev.net/
   4
   5 The Configurable Math Library (CML) is distributed under the terms of the
   6 Boost Software License, v1.0 (see cml/LICENSE for details).
   7
   8  *-----------------------------------------------------------------------*/
   9 /** @file
  10  *  @brief
  11  */
  12
  13 #ifndef coord_conversion_h
  14 #define coord_conversion_h
  15
  16 #include <cml/mathlib/checking.h>
  17 #include <cml/mathlib/epsilon.h>
  18 #include <cml/mathlib/helper.h>
  19
  20 /* Functions for converting between Cartesian, polar, cylindrical and
  21  * spherical coordinates.
  22  *
  23  * The 3D conversion functions take an integer axis index argument. For
  24  * cylindrical coordinates this determines the axis of the cylinder, and for
  25  * spherical it determines which cardinal axis is normal to the azimuth plane.
  26  *
  27  * For spherical coordinates the option of whether to treat phi as latitude
  28  * or colatitude is also available. The 'type' argument takes either of the
  29  * enumerants cml::latitude and cml::colatitude to reflect this.
  30  */
  31
  32 namespace cml {
  33
  34 //////////////////////////////////////////////////////////////////////////////
  35 // Conversion to Cartesian coordinates
  36 //////////////////////////////////////////////////////////////////////////////
  37
  38 /* Convert cylindrical coordinates to Cartesian coordinates in R3 */
  39 template < typename E, class A > void
  40 cylindrical_to_cartesian(
  41     E radius, E theta, E height, size_t axis, vector<E,A>& v)
  42 {
  43     typedef vector<E,A> vector_type;
  44     typedef typename vector_type::value_type value_type;
  45
  46     /* Checking */
  47     detail::CheckVec3(v);
  48     detail::CheckIndex3(axis);
  49
  50     size_t i, j, k;
  51     cyclic_permutation(axis, i, j, k);
  52
  53     v[i] = height;
  54     v[j] = std::cos(theta) * radius;
  55     v[k] = std::sin(theta) * radius;
  56 }
  57
  58 /* Convert spherical coordinates to Cartesian coordinates in R3 */
  59 template < typename E, class A > void
  60 spherical_to_cartesian(E radius, E theta, E phi, size_t axis,
  61     SphericalType type, vector<E,A>& v)
  62 {
  63     typedef vector<E,A> vector_type;
  64     typedef typename vector_type::value_type value_type;
  65
  66     /* Checking */
  67     detail::CheckVec3(v);
  68     detail::CheckIndex3(axis);
  69
  70     if (type == latitude) {
  71         phi = constants<value_type>::pi_over_2() - phi;
  72     }
  73
  74     value_type sin_phi = std::sin(phi);
  75     value_type cos_phi = std::cos(phi);
  76     value_type sin_phi_r = sin_phi * radius;
  77
  78     size_t i, j, k;
  79     cyclic_permutation(axis, i, j, k);
  80
  81     v[i] = cos_phi * radius;
  82     v[j] = sin_phi_r * std::cos(theta);
  83     v[k] = sin_phi_r * std::sin(theta);
  84 }
  85
  86 /* Convert polar coordinates to Cartesian coordinates in R2 */
  87 template < typename E, class A > void
  88 polar_to_cartesian(E radius, E theta, vector<E,A>& v)
  89 {
  90     /* Checking handled by set() */
  91     v.set(std::cos(theta) * double(radius), std::sin(theta) * double(radius));
  92 }
  93
  94 //////////////////////////////////////////////////////////////////////////////
  95 // Conversion from Cartesian coordinates
  96 //////////////////////////////////////////////////////////////////////////////
  97
  98 /* Convert Cartesian coordinates to cylindrical coordinates in R3  */
  99 template < class VecT, typename Real > void
 100 cartesian_to_cylindrical(const VecT& v, Real& radius, Real& theta,
 101     Real& height, size_t axis, Real tolerance = epsilon<Real>::placeholder())
 102 {
 103     typedef Real value_type;
 104
 105     /* Checking */
 106     detail::CheckVec3(v);
 107     detail::CheckIndex3(axis);
 108
 109     size_t i, j, k;
 110     cyclic_permutation(axis, i, j, k);
 111
 112     radius = length(v[j],v[k]);
 113     theta = radius < tolerance ? value_type(0) : std::atan2(v[k],v[j]);
 114     height = v[i];
 115 }
 116
 117 /* Convert Cartesian coordinates to spherical coordinates in R3 */
 118 template < class VecT, typename Real > void
 119 cartesian_to_spherical(const VecT& v, Real& radius, Real& theta, Real& phi,
 120     size_t axis, SphericalType type,
 121     Real tolerance = epsilon<Real>::placeholder())
 122 {
 123     typedef Real value_type;
 124
 125     /* Checking */
 126     detail::CheckVec3(v);
 127     detail::CheckIndex3(axis);
 128
 129     size_t i, j, k;
 130     cyclic_permutation(axis, i, j, k);
 131
 132     value_type len = length(v[j],v[k]);
 133     theta = len < tolerance ? value_type(0) : std::atan2(v[k],v[j]);
 134     radius = length(v[i], len);
 135     if (radius < tolerance) {
 136         phi = value_type(0);
 137     } else {
 138         phi = std::atan2(len,v[i]);
 139         //phi = type.convert(phi);
 140         if (type == latitude) {
 141             phi = constants<value_type>::pi_over_2() - phi;
 142         }
 143     }
 144 }
 145
 146 /* Convert Cartesian coordinates to polar coordinates in R2 */
 147 template < class VecT, typename Real > void
 148 cartesian_to_polar(const VecT& v, Real& radius, Real& theta,
 149     Real tolerance = epsilon<Real>::placeholder())
 150 {
 151     typedef Real value_type;
 152
 153     /* Checking */
 154     detail::CheckVec2(v);
 155
 156     radius = v.length();
 157     theta = radius < tolerance ? value_type(0) : std::atan2(v[1],v[0]);
 158 }
 159
 160 } // namespace cml
 161
 162 #endif