Line data Source code
1 : /** \file zernike.cpp
2 : * \author Jared R. Males (jaredmales@gmail.com)
3 : * \brief Working with the Zernike polynomials.
4 : *
5 : * \ingroup signal_processing_files
6 : *
7 : */
8 :
9 : //***********************************************************************//
10 : // Copyright 2021 Jared R. Males (jaredmales@gmail.com)
11 : //
12 : // This file is part of mxlib.
13 : //
14 : // mxlib is free software: you can redistribute it and/or modify
15 : // it under the terms of the GNU General Public License as published by
16 : // the Free Software Foundation, either version 3 of the License, or
17 : // (at your option) any later version.
18 : //
19 : // mxlib is distributed in the hope that it will be useful,
20 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
21 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 : // GNU General Public License for more details.
23 : //
24 : // You should have received a copy of the GNU General Public License
25 : // along with mxlib. If not, see <http://www.gnu.org/licenses/>.
26 : //***********************************************************************//
27 :
28 : #include "sigproc/zernike.hpp"
29 :
30 : namespace mx
31 : {
32 : namespace sigproc
33 : {
34 :
35 28 : int noll_nm( int &n, int &m, int j )
36 : {
37 28 : if( j < 1 )
38 : {
39 2 : internal::mxlib_error_report(error_t::invalidarg, "The Noll index j cannot be less than 1 in the Zernike polynomials" );
40 1 : return -1;
41 : }
42 :
43 27 : n = ceil( -1.5 + sqrt( 0.25 + 2 * j ) - 1e-10 ); // 1e-10 is to avoid wrong rounding due to numerical precision
44 :
45 27 : int jrem = j - ( n * ( n + 1 ) / 2 + 1 );
46 54 : m = ( jrem + ( jrem % 2 ) * abs( ( n % 2 ) - 1 ) + fabs( ( jrem % 2 ) - 1 ) * ( n % 2 ) ) *
47 27 : ( -math::func::sign( ( j % 2 ) - 0.5 ) );
48 :
49 27 : return 0;
50 : }
51 :
52 0 : int noll_j( unsigned int n, int m )
53 : {
54 0 : if( ( ( n - m ) % 2 ) )
55 0 : return -1; // Check if odd
56 :
57 0 : int mn = n % 4;
58 :
59 0 : int dm = 0;
60 0 : if( m >= 0 && ( mn == 2 || mn == 3 ) )
61 0 : dm = 1;
62 0 : else if( m <= 0 && ( mn == 0 || mn == 1 ) )
63 0 : dm = 1;
64 :
65 0 : int j = ( n * ( n + 1 ) ) / 2 + abs( m ) + dm;
66 :
67 0 : return j;
68 : }
69 :
70 0 : int nZernRadOrd( unsigned int n )
71 : {
72 0 : if( n % 2 ) // odd
73 : {
74 0 : return noll_j( n, -n );
75 : }
76 : else
77 : {
78 0 : return noll_j( n, n );
79 : }
80 : }
81 :
82 : // Explicit instantiations:
83 : template int zernikeRCoeffs<float>( std::vector<float> &c, int n, int m );
84 :
85 : template int zernikeRCoeffs<double>( std::vector<double> &c, int n, int m );
86 :
87 : template int zernikeRCoeffs<long double>( std::vector<long double> &c, int n, int m );
88 :
89 : #ifdef HASQUAD
90 : template int zernikeRCoeffs<__float128>( std::vector<__float128> &c, int n, int m );
91 : #endif
92 :
93 : template float zernikeR<float, double>( float rho, int n, int m, std::vector<double> &c );
94 :
95 : template double zernikeR<double, double>( double rho, int n, int m, std::vector<double> &c );
96 :
97 : template long double zernikeR<long double, long double>( long double rho, int n, int m, std::vector<long double> &c );
98 :
99 : #ifdef HASQUAD
100 : template __float128 zernikeR<__float128, __float128>( __float128 rho, int n, int m, std::vector<__float128> &c );
101 : #endif
102 :
103 : template float zernikeR<float, double>( float rho, int n, int m );
104 :
105 : template double zernikeR<double, double>( double rho, int n, int m );
106 :
107 : template long double zernikeR<long double, long double>( long double rho, int n, int m );
108 :
109 : #ifdef HASQUAD
110 : template __float128 zernikeR<__float128, __float128>( __float128 rho, int n, int m );
111 : #endif
112 :
113 : template float zernike<float, double>( float rho, float phi, int n, int m, std::vector<double> &c );
114 :
115 : template double zernike<double, double>( double rho, double phi, int n, int m, std::vector<double> &c );
116 :
117 : template long double
118 : zernike<long double, long double>( long double rho, long double phi, int n, int m, std::vector<long double> &c );
119 :
120 : #ifdef HASQUAD
121 : template __float128
122 : zernike<__float128, __float128>( __float128 rho, __float128 phi, int n, int m, std::vector<__float128> &c );
123 : #endif
124 :
125 : template float zernike<float, double>( float rho, float phi, int n, int m );
126 :
127 : template double zernike<double, double>( double rho, double phi, int n, int m );
128 :
129 : template long double zernike<long double, long double>( long double rho, long double phi, int n, int m );
130 :
131 : #ifdef HASQUAD
132 : template __float128 zernike<__float128, __float128>( __float128 rho, __float128 phi, int n, int m );
133 : #endif
134 :
135 : template float zernike<float, double>( float rho, float phi, int j );
136 :
137 : template double zernike<double, double>( double rho, double phi, int j );
138 :
139 : template long double zernike<long double, long double>( long double rho, long double phi, int j );
140 :
141 : #ifdef HASQUAD
142 : template __float128 zernike<__float128, __float128>( __float128 rho, __float128 phi, int j );
143 : #endif
144 :
145 : template float zernikeQNorm<float>( float k, float phi, int n, int m );
146 :
147 : template double zernikeQNorm<double>( double k, double phi, int n, int m );
148 :
149 : template long double zernikeQNorm<long double>( long double k, long double phi, int n, int m );
150 :
151 : #ifdef HASQUAD
152 : template __float128 zernikeQNorm<__float128>( __float128 k, __float128 phi, int n, int m );
153 : #endif
154 :
155 : } // namespace sigproc
156 : } // namespace mx
|