Planet
navi homePPSaboutscreenshotsdownloaddevelopmentforum

source: downloads/boost_1_34_1/libs/math/quaternion/HSO3SO4.cpp @ 33

Last change on this file since 33 was 29, checked in by landauf, 16 years ago

updated boost from 1_33_1 to 1_34_1

File size: 16.1 KB
Line 
1// test file for HSO3.hpp and HSO4.hpp
2
3//  (C) Copyright Hubert Holin 2001.
4//  Distributed under the Boost Software License, Version 1.0. (See
5//  accompanying file LICENSE_1_0.txt or copy at
6//  http://www.boost.org/LICENSE_1_0.txt)
7
8
9#include <iostream>
10
11#include <boost/math/quaternion.hpp>
12
13#include "HSO3.hpp"
14#include "HSO4.hpp"
15
16
17const int    number_of_intervals = 5;
18
19const float    pi = ::std::atan(1.0f)*4;
20
21
22
23void    test_SO3();
24   
25void    test_SO4();
26
27
28int    main()
29
30{
31    test_SO3();
32   
33    test_SO4();
34   
35    ::std::cout << "That's all folks!" << ::std::endl;
36}
37
38
39//
40//    Test of quaternion and R^3 rotation relationship
41//
42
43void    test_SO3_spherical()
44{
45    ::std::cout << "Testing spherical:" << ::std::endl;
46    ::std::cout << ::std::endl;
47   
48    const float    rho = 1.0f;
49   
50    float        theta;
51    float        phi1;
52    float        phi2;
53   
54    for    (int idxphi2 = 0; idxphi2 <= number_of_intervals; idxphi2++)
55    {
56        phi2 = (-pi/2)+(idxphi2*pi)/number_of_intervals;
57       
58        for    (int idxphi1 = 0; idxphi1 <= number_of_intervals; idxphi1++)
59        {
60            phi1 = (-pi/2)+(idxphi1*pi)/number_of_intervals;
61           
62            for    (int idxtheta = 0; idxtheta <= number_of_intervals; idxtheta++)
63            {
64                theta = -pi+(idxtheta*(2*pi))/number_of_intervals;
65               
66                //::std::cout << "theta = " << theta << " ; ";
67                //::std::cout << "phi1 = " << phi1 << " ; ";
68                //::std::cout << "phi2 = " << phi2;
69                //::std::cout << ::std::endl;
70               
71                ::boost::math::quaternion<float>    q = ::boost::math::spherical(rho, theta, phi1, phi2);
72               
73                //::std::cout << "q = " << q << ::std::endl;
74               
75                R3_matrix<float>                    rot = quaternion_to_R3_rotation(q);
76               
77                //::std::cout << "rot = ";
78                //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << ::std::endl;
79                //::std::cout << "\t";
80                //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << ::std::endl;
81                //::std::cout << "\t";
82                //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << ::std::endl;
83               
84                ::boost::math::quaternion<float>    p = R3_rotation_to_quaternion(rot, &q);
85               
86                //::std::cout << "p = " << p << ::std::endl;
87               
88                //::std::cout << "round trip discrepancy: " << ::boost::math::abs(q-p) << ::std::endl;
89               
90                //::std::cout << ::std::endl;
91            }
92        }
93    }
94   
95    ::std::cout << ::std::endl;
96}
97
98   
99void    test_SO3_semipolar()
100{
101    ::std::cout << "Testing semipolar:" << ::std::endl;
102    ::std::cout << ::std::endl;
103   
104    const float    rho = 1.0f;
105   
106    float        alpha;
107    float        theta1;
108    float        theta2;
109   
110    for    (int idxalpha = 0; idxalpha <= number_of_intervals; idxalpha++)
111    {
112        alpha = (idxalpha*(pi/2))/number_of_intervals;
113       
114        for    (int idxtheta1 = 0; idxtheta1 <= number_of_intervals; idxtheta1++)
115        {
116            theta1 = -pi+(idxtheta1*(2*pi))/number_of_intervals;
117           
118            for    (int idxtheta2 = 0; idxtheta2 <= number_of_intervals; idxtheta2++)
119            {
120                theta2 = -pi+(idxtheta2*(2*pi))/number_of_intervals;
121               
122                //::std::cout << "alpha = " << alpha << " ; ";
123                //::std::cout << "theta1 = " << theta1 << " ; ";
124                //::std::cout << "theta2 = " << theta2;
125                //::std::cout << ::std::endl;
126               
127                ::boost::math::quaternion<float>    q = ::boost::math::semipolar(rho, alpha, theta1, theta2);
128               
129                //::std::cout << "q = " << q << ::std::endl;
130               
131                R3_matrix<float>                    rot = quaternion_to_R3_rotation(q);
132               
133                //::std::cout << "rot = ";
134                //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << ::std::endl;
135                //::std::cout << "\t";
136                //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << ::std::endl;
137                //::std::cout << "\t";
138                //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << ::std::endl;
139               
140                ::boost::math::quaternion<float>    p = R3_rotation_to_quaternion(rot, &q);
141               
142                //::std::cout << "p = " << p << ::std::endl;
143               
144                //::std::cout << "round trip discrepancy: " << ::boost::math::abs(q-p) << ::std::endl;
145               
146                //::std::cout << ::std::endl;
147            }
148        }
149    }
150   
151    ::std::cout << ::std::endl;
152}
153
154   
155void    test_SO3_multipolar()
156{
157    ::std::cout << "Testing multipolar:" << ::std::endl;
158    ::std::cout << ::std::endl;
159   
160    float    rho1;
161    float    rho2;
162   
163    float    theta1;
164    float    theta2;
165   
166    for    (int idxrho = 0; idxrho <= number_of_intervals; idxrho++)
167    {
168        rho1 = (idxrho*1.0f)/number_of_intervals;
169        rho2 = ::std::sqrt(1.0f-rho1*rho1);
170       
171        for    (int idxtheta1 = 0; idxtheta1 <= number_of_intervals; idxtheta1++)
172        {
173            theta1 = -pi+(idxtheta1*(2*pi))/number_of_intervals;
174           
175            for    (int idxtheta2 = 0; idxtheta2 <= number_of_intervals; idxtheta2++)
176            {
177                theta2 = -pi+(idxtheta2*(2*pi))/number_of_intervals;
178               
179                //::std::cout << "rho1 = " << rho1 << " ; ";
180                //::std::cout << "theta1 = " << theta1 << " ; ";
181                //::std::cout << "theta2 = " << theta2;
182                //::std::cout << ::std::endl;
183               
184                ::boost::math::quaternion<float>    q = ::boost::math::multipolar(rho1, theta1, rho2, theta2);
185               
186                //::std::cout << "q = " << q << ::std::endl;
187               
188                R3_matrix<float>                    rot = quaternion_to_R3_rotation(q);
189               
190                //::std::cout << "rot = ";
191                //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << ::std::endl;
192                //::std::cout << "\t";
193                //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << ::std::endl;
194                //::std::cout << "\t";
195                //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << ::std::endl;
196               
197                ::boost::math::quaternion<float>    p = R3_rotation_to_quaternion(rot, &q);
198               
199                //::std::cout << "p = " << p << ::std::endl;
200               
201                //::std::cout << "round trip discrepancy: " << ::boost::math::abs(q-p) << ::std::endl;
202               
203                //::std::cout << ::std::endl;
204            }
205        }
206    }
207   
208    ::std::cout << ::std::endl;
209}
210
211   
212void    test_SO3_cylindrospherical()
213{
214    ::std::cout << "Testing cylindrospherical:" << ::std::endl;
215    ::std::cout << ::std::endl;
216   
217    float    t;
218   
219    float    radius;
220    float    longitude;
221    float    latitude;
222   
223    for    (int idxt = 0; idxt <= number_of_intervals; idxt++)
224    {
225        t = -1.0f+(idxt*2.0f)/number_of_intervals;
226        radius = ::std::sqrt(1.0f-t*t);
227       
228        for    (int idxlatitude = 0; idxlatitude <= number_of_intervals; idxlatitude++)
229        {
230            latitude = (-pi/2)+(idxlatitude*pi)/number_of_intervals;
231           
232            for    (int idxlongitude = 0; idxlongitude <= number_of_intervals; idxlongitude++)
233            {
234                longitude = -pi+(idxlongitude*(2*pi))/number_of_intervals;
235               
236                //::std::cout << "t = " << t << " ; ";
237                //::std::cout << "longitude = " << longitude;
238                //::std::cout << "latitude = " << latitude;
239                //::std::cout << ::std::endl;
240               
241                ::boost::math::quaternion<float>    q = ::boost::math::cylindrospherical(t, radius, longitude, latitude);
242               
243                //::std::cout << "q = " << q << ::std::endl;
244               
245                R3_matrix<float>                    rot = quaternion_to_R3_rotation(q);
246               
247                //::std::cout << "rot = ";
248                //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << ::std::endl;
249                //::std::cout << "\t";
250                //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << ::std::endl;
251                //::std::cout << "\t";
252                //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << ::std::endl;
253               
254                ::boost::math::quaternion<float>    p = R3_rotation_to_quaternion(rot, &q);
255               
256                //::std::cout << "p = " << p << ::std::endl;
257               
258                //::std::cout << "round trip discrepancy: " << ::boost::math::abs(q-p) << ::std::endl;
259               
260                //::std::cout << ::std::endl;
261            }
262        }
263    }
264   
265    ::std::cout << ::std::endl;
266}
267
268   
269void    test_SO3_cylindrical()
270{
271    ::std::cout << "Testing cylindrical:" << ::std::endl;
272    ::std::cout << ::std::endl;
273   
274    float    r;
275    float    angle;
276   
277    float    h1;
278    float    h2;
279   
280    for    (int idxh2 = 0; idxh2 <= number_of_intervals; idxh2++)
281    {
282        h2 = -1.0f+(idxh2*2.0f)/number_of_intervals;
283       
284        for    (int idxh1 = 0; idxh1 <= number_of_intervals; idxh1++)
285        {
286            h1 = ::std::sqrt(1.0f-h2*h2)*(-1.0f+(idxh2*2.0f)/number_of_intervals);
287            r = ::std::sqrt(1.0f-h1*h1-h2*h2);
288           
289            for    (int idxangle = 0; idxangle <= number_of_intervals; idxangle++)
290            {
291                angle = -pi+(idxangle*(2*pi))/number_of_intervals;
292               
293                //::std::cout << "angle = " << angle << " ; ";
294                //::std::cout << "h1 = " << h1;
295                //::std::cout << "h2 = " << h2;
296                //::std::cout << ::std::endl;
297               
298                ::boost::math::quaternion<float>    q = ::boost::math::cylindrical(r, angle, h1, h2);
299               
300                //::std::cout << "q = " << q << ::std::endl;
301               
302                R3_matrix<float>                    rot = quaternion_to_R3_rotation(q);
303               
304                //::std::cout << "rot = ";
305                //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << ::std::endl;
306                //::std::cout << "\t";
307                //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << ::std::endl;
308                //::std::cout << "\t";
309                //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << ::std::endl;
310               
311                ::boost::math::quaternion<float>    p = R3_rotation_to_quaternion(rot, &q);
312               
313                //::std::cout << "p = " << p << ::std::endl;
314               
315                //::std::cout << "round trip discrepancy: " << ::boost::math::abs(q-p) << ::std::endl;
316               
317                //::std::cout << ::std::endl;
318            }
319        }
320    }
321   
322    ::std::cout << ::std::endl;
323}
324
325
326void    test_SO3()
327{
328    ::std::cout << "Testing SO3:" << ::std::endl;
329    ::std::cout << ::std::endl;
330   
331    test_SO3_spherical();
332   
333    test_SO3_semipolar();
334   
335    test_SO3_multipolar();
336   
337    test_SO3_cylindrospherical();
338   
339    test_SO3_cylindrical();
340}
341
342
343//
344//    Test of quaternion and R^4 rotation relationship
345//
346
347void    test_SO4_spherical()
348{
349    ::std::cout << "Testing spherical:" << ::std::endl;
350    ::std::cout << ::std::endl;
351   
352    const float    rho1 = 1.0f;
353    const float    rho2 = 1.0f;
354   
355    float        theta1;
356    float        phi11;
357    float        phi21;
358   
359    float        theta2;
360    float        phi12;
361    float        phi22;
362   
363    for    (int idxphi21 = 0; idxphi21 <= number_of_intervals; idxphi21++)
364    {
365        phi21 = (-pi/2)+(idxphi21*pi)/number_of_intervals;
366       
367        for    (int idxphi22 = 0; idxphi22 <= number_of_intervals; idxphi22++)
368        {
369            phi22 = (-pi/2)+(idxphi22*pi)/number_of_intervals;
370           
371            for    (int idxphi11 = 0; idxphi11 <= number_of_intervals; idxphi11++)
372            {
373                phi11 = (-pi/2)+(idxphi11*pi)/number_of_intervals;
374               
375                for    (int idxphi12 = 0; idxphi12 <= number_of_intervals; idxphi12++)
376                {
377                    phi12 = (-pi/2)+(idxphi12*pi)/number_of_intervals;
378                   
379                    for    (int idxtheta1 = 0; idxtheta1 <= number_of_intervals; idxtheta1++)
380                    {
381                        theta1 = -pi+(idxtheta1*(2*pi))/number_of_intervals;
382                       
383                        for    (int idxtheta2 = 0; idxtheta2 <= number_of_intervals; idxtheta2++)
384                        {
385                            theta2 = -pi+(idxtheta2*(2*pi))/number_of_intervals;
386                           
387                            //::std::cout << "theta1 = " << theta1 << " ; ";
388                            //::std::cout << "phi11 = " << phi11 << " ; ";
389                            //::std::cout << "phi21 = " << phi21;
390                            //::std::cout << "theta2 = " << theta2 << " ; ";
391                            //::std::cout << "phi12 = " << phi12 << " ; ";
392                            //::std::cout << "phi22 = " << phi22;
393                            //::std::cout << ::std::endl;
394                           
395                            ::boost::math::quaternion<float>    p1 = ::boost::math::spherical(rho1, theta1, phi11, phi21);
396                           
397                            //::std::cout << "p1 = " << p1 << ::std::endl;
398                           
399                            ::boost::math::quaternion<float>    q1 = ::boost::math::spherical(rho2, theta2, phi12, phi22);
400                           
401                            //::std::cout << "q1 = " << q1 << ::std::endl;
402                           
403                            ::std::pair< ::boost::math::quaternion<float> , ::boost::math::quaternion<float> >    pq1 =
404                                ::std::make_pair(p1,q1);
405                           
406                            R4_matrix<float>                    rot = quaternions_to_R4_rotation(pq1);
407                           
408                            //::std::cout << "rot = ";
409                            //::std::cout << "\t" << rot.a11 << "\t" << rot.a12 << "\t" << rot.a13 << "\t" << rot.a14 << ::std::endl;
410                            //::std::cout << "\t";
411                            //::std::cout << "\t" << rot.a21 << "\t" << rot.a22 << "\t" << rot.a23 << "\t" << rot.a24 << ::std::endl;
412                            //::std::cout << "\t";
413                            //::std::cout << "\t" << rot.a31 << "\t" << rot.a32 << "\t" << rot.a33 << "\t" << rot.a34 << ::std::endl;
414                            //::std::cout << "\t";
415                            //::std::cout << "\t" << rot.a41 << "\t" << rot.a42 << "\t" << rot.a43 << "\t" << rot.a44 << ::std::endl;
416                           
417                            ::std::pair< ::boost::math::quaternion<float> , ::boost::math::quaternion<float> >    pq2 =
418                                R4_rotation_to_quaternions(rot, &pq1);
419                           
420                            //::std::cout << "p1 = " << pq.first << ::std::endl;
421                            //::std::cout << "p2 = " << pq.second << ::std::endl;
422                           
423                            //::std::cout << "round trip discrepancy: " << ::std::sqrt(::boost::math::norm(pq1.first-pq2.first)+::boost::math::norm(pq1.second-pq2.second)) << ::std::endl;
424                           
425                            //::std::cout << ::std::endl;
426                        }
427                    }
428                }
429            }
430        }
431    }
432   
433    ::std::cout << ::std::endl;
434}
435
436
437void    test_SO4()
438{
439    ::std::cout << "Testing SO4:" << ::std::endl;
440    ::std::cout << ::std::endl;
441   
442    test_SO4_spherical();
443}
444
445
Note: See TracBrowser for help on using the repository browser.