1 /*
2 * Copyright (c) 2011-2014, The University of Oxford
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 * 1. Redistributions of source code must retain the above copyright notice,
8 * this list of conditions and the following disclaimer.
9 * 2. Redistributions in binary form must reproduce the above copyright notice,
10 * this list of conditions and the following disclaimer in the documentation
11 * and/or other materials provided with the distribution.
12 * 3. Neither the name of the University of Oxford nor the names of its
13 * contributors may be used to endorse or promote products derived from this
14 * software without specific prior written permission.
15 *
16 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
19 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
20 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
21 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
22 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
23 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
24 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
25 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
26 * POSSIBILITY OF SUCH DAMAGE.
27 */
28
29
30 #include "math/oskar_sph_rotate_to_position.h"
31 #include <stdlib.h>
32 #include <math.h>
33
34 #ifdef __cplusplus
35 extern "C" {
36 #endif
37
oskar_sph_rotate_to_position(int n,oskar_Mem * lon,oskar_Mem * lat,double lon0,double lat0)38 int oskar_sph_rotate_to_position(int n, oskar_Mem* lon, oskar_Mem* lat,
39 double lon0, double lat0)
40 {
41 int i;
42
43 if (lon == NULL || lat == NULL)
44 return OSKAR_ERR_INVALID_ARGUMENT;
45
46 if (oskar_mem_location(lon) != OSKAR_CPU || oskar_mem_location(lat) != OSKAR_CPU)
47 return OSKAR_ERR_BAD_LOCATION;
48
49 if (oskar_mem_type(lon) == OSKAR_DOUBLE && oskar_mem_type(lat) == OSKAR_DOUBLE)
50 {
51 /* Construct rotation matrix. */
52 double cosLon, sinLon, cosLat, sinLat;
53 double m11, m12, m13, m21, m22, m23, m31, m33;
54 cosLon = cos(lon0);
55 sinLon = sin(lon0);
56 cosLat = cos(lat0);
57 sinLat = sin(lat0);
58 m11 = cosLon * sinLat; m12 = -sinLon; m13 = cosLon * cosLat;
59 m21 = sinLon * sinLat; m22 = cosLon; m23 = sinLon * cosLat;
60 m31 = -cosLat; m33 = sinLat;
61
62 for (i = 0; i < n; ++i)
63 {
64 double x, y, z, a, b, c;
65
66 /* Direction cosines */
67 c = cos(((double*)lat->data)[i]);
68 x = c * cos(((double*)lon->data)[i]);
69 y = c * sin(((double*)lon->data)[i]);
70 z = sin(((double*)lat->data)[i]);
71
72 /* Apply rotation matrix */
73 a = m11 * x + m12 * y + m13 * z;
74 b = m21 * x + m22 * y + m23 * z;
75 c = m31 * x + m33 * z;
76
77 /* Convert back to angles. */
78 ((double*)lon->data)[i] = atan2(b, a);
79 ((double*)lat->data)[i] = atan2(c, sqrt(a*a + b*b));
80 }
81 }
82 else if (oskar_mem_type(lon) == OSKAR_SINGLE && oskar_mem_type(lat) == OSKAR_SINGLE)
83 {
84 /* Construct rotation matrix. */
85 float cosLon, sinLon, cosLat, sinLat;
86 float m11, m12, m13, m21, m22, m23, m31, m33;
87 cosLon = cosf(lon0);
88 sinLon = sinf(lon0);
89 cosLat = cosf(lat0);
90 sinLat = sinf(lat0);
91 m11 = cosLon * sinLat; m12 = -sinLon; m13 = cosLon * cosLat;
92 m21 = sinLon * sinLat; m22 = cosLon; m23 = sinLon * cosLat;
93 m31 = -cosLat; m33 = sinLat;
94
95 for (i = 0; i < n; ++i)
96 {
97 float x, y, z, a, b, c;
98
99 /* Direction cosines */
100 c = cosf(((float*)lat->data)[i]);
101 x = c * cosf(((float*)lon->data)[i]);
102 y = c * sinf(((float*)lon->data)[i]);
103 z = sinf(((float*)lat->data)[i]);
104
105 /* Apply rotation matrix */
106 a = m11 * x + m12 * y + m13 * z;
107 b = m21 * x + m22 * y + m23 * z;
108 c = m31 * x + m33 * z;
109
110 /* Convert back to angles. */
111 ((float*)lon->data)[i] = atan2f(b, a);
112 ((float*)lat->data)[i] = atan2f(c, sqrtf(a*a + b*b));
113 }
114 }
115 else
116 {
117 return OSKAR_ERR_BAD_DATA_TYPE;
118 }
119
120 return 0;
121 }
122
123
124 #ifdef __cplusplus
125 }
126 #endif
127