~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

Linux Cross Reference
Tina5/tina-libs/tina/math/mathTyp_cmplx.c

Version: ~
Architecture: ~ [ i386 ] ~ [ alpha ] ~ [ m68k ] ~ [ mips ] ~ [ ppc ] ~ [ sparc ] ~ [ sparc64 ] ~

  1 /**********
  2  *
  3  * Copyright (c) 2003, Division of Imaging Science and Biomedical Engineering,
  4  * University of Manchester, UK.  All rights reserved.
  5  * 
  6  * Redistribution and use in source and binary forms, with or without modification, 
  7  * are permitted provided that the following conditions are met:
  8  * 
  9  *   . Redistributions of source code must retain the above copyright notice, 
 10  *     this list of conditions and the following disclaimer.
 11  *    
 12  *   . Redistributions in binary form must reproduce the above copyright notice,
 13  *     this list of conditions and the following disclaimer in the documentation 
 14  *     and/or other materials provided with the distribution.
 15  * 
 16  *   . Neither the name of the University of Manchester nor the names of its
 17  *     contributors may be used to endorse or promote products derived from this 
 18  *     software without specific prior written permission.
 19  * 
 20  * 
 21  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" 
 22  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 
 23  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 
 24  * ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE 
 25  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
 26  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
 27  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
 28  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 
 29  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) 
 30  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 
 31  * POSSIBILITY OF SUCH DAMAGE.
 32  *
 33  **********
 34  *
 35  * Program :    TINA
 36  * File    :  $Source: /home/tina/cvs/tina-libs/tina/math/mathTyp_cmplx.c,v $
 37  * Date    :  $Date: 2003/09/23 11:32:20 $
 38  * Version :  $Revision: 1.4 $
 39  * CVS Id  :  $Id: mathTyp_cmplx.c,v 1.4 2003/09/23 11:32:20 matts Exp $
 40  *
 41  * Notes :
 42  *
 43  *  Complex handling.
 44  *  Complex is {float x; float y;}
 45  *
 46  *********
 47 */
 48 /** 
 49  *  @file
 50  *  @brief  Complex arithmetic functions.       
 51  *
 52  *  Complex is {float x; float y;}.  
 53  * 
 54 */
 55 
 56 #include "mathTyp_cmplx.h"
 57 
 58 #if HAVE_CONFIG_H
 59 #include <config.h>
 60 #endif
 61 
 62 #include <math.h>
 63 
 64 Complex         cmplx_sum(Complex a, Complex b)
 65 {
 66         Complex         c = {Complex_id};
 67 
 68         c.x = a.x + b.x;
 69         c.y = a.y + b.y;
 70         return (c);
 71 }
 72 
 73 Complex         cmplx_diff(Complex a, Complex b)
 74 {
 75         Complex         c = {Complex_id};
 76 
 77         c.x = a.x - b.x;
 78         c.y = a.y - b.y;
 79         return (c);
 80 }
 81 
 82 Complex         cmplx_minus(Complex a)
 83 {
 84         Complex         c = {Complex_id};
 85 
 86         c.x = -a.x;
 87         c.y = -a.y;
 88         return (c);
 89 }
 90 
 91 Complex         cmplx_prod(Complex a, Complex b)
 92 {
 93         Complex         c = {Complex_id};
 94 
 95         c.x = a.x * b.x - a.y * b.y;
 96         c.y = a.y * b.x + a.x * b.y;
 97         return (c);
 98 }
 99 
100 Complex         cmplx_cprod(Complex a, Complex b)
101 {
102         Complex         c = {Complex_id};
103 
104         /** a times conjugate of b **/
105         c.x = a.x * b.x + a.y * b.y;
106         c.y = a.y * b.x - a.x * b.y;
107         return (c);
108 }
109 
110 Complex         cmplx_times(double a, Complex b)
111 {
112         Complex         c = {Complex_id};
113 
114         c.x = a * b.x;
115         c.y = a * b.y;
116         return (c);
117 }
118 
119 Complex         cmplx(double x, double y)
120 {
121         Complex         c = {Complex_id};
122 
123         c.x = x;
124         c.y = y;
125         return (c);
126 }
127 
128 Complex         cmplx_zero(void)
129 {
130         Complex         zero = {Complex_id};
131 
132         zero.x = zero.y = 0.0;
133         return (zero);
134 }
135 
136 Complex         cmplx_unit(void)
137 {
138         Complex         unit = {Complex_id};
139 
140         unit.x = 1.0;
141         unit.y = 0.0;
142         return (unit);
143 }
144 
145 Complex         cmplx_conj(Complex z)
146 {
147         Complex         c = {Complex_id};
148 
149         c.x = z.x;
150         c.y = -z.y;
151         return (c);
152 }
153 
154 Complex         cmplx_div(Complex a, Complex b)
155 {
156         Complex         c = {Complex_id};
157         double          r, den;
158 
159         if (fabs(b.x) >= fabs(b.y))
160         {
161                 r = b.y / b.x;
162                 den = b.x + r * b.y;
163                 c.x = (a.x + r * a.y) / den;
164                 c.y = (a.y - r * a.x) / den;
165         } else
166         {
167                 r = b.x / b.y;
168                 den = b.y + r * b.x;
169                 c.x = (a.x * r + a.y) / den;
170                 c.y = (a.y * r - a.x) / den;
171         }
172         return (c);
173 }
174 
175 double          cmplx_mod(Complex z)
176 {
177         double          x, y, ans, temp;
178 
179         x = fabs(z.x);
180         y = fabs(z.y);
181         if (x == 0.0)
182                 ans = y;
183         else if (y == 0.0)
184                 ans = x;
185         else if (x > y)
186         {
187                 temp = y / x;
188                 ans = x * sqrt(1.0 + temp * temp);
189         } else
190         {
191                 temp = x / y;
192                 ans = y * sqrt(1.0 + temp * temp);
193         }
194         return (ans);
195 }
196 
197 double          cmplx_sqrmod(Complex z)
198 {
199         double          x, y;
200 
201         x = z.x;
202         y = z.y;
203         return (x * x + y * y);
204 }
205 
206 double          cmplx_arg(Complex z)
207 {
208         if (z.x == 0.0 && z.y == 0.0)
209                 return (0.0);
210         return (atan2(z.y, z.x));
211 }
212 
213 Complex         cmplx_phase(Complex z)
214 {
215         double          r = cmplx_mod(z);
216 
217         if (r == 0.0)
218                 return (cmplx(1.0, 0.0));
219         else
220                 return (cmplx(z.x / r, z.y / r));
221 }
222 
223 Complex         cmplx_sqr(Complex a)
224 {
225         Complex         c = {Complex_id};
226 
227         c.x = a.x * a.x - a.y * a.y;
228         c.y = a.y * a.x + a.x * a.y;
229         return (c);
230 }
231 
232 Complex         cmplx_sqrt(Complex z)
233 {
234         Complex         c = {Complex_id};
235         double          x, y, w, r;
236 
237         if ((z.x == 0.0) && (z.y == 0.0))
238         {
239                 c.x = 0.0;
240                 c.y = 0.0;
241                 return (c);
242         } else
243         {
244                 x = fabs(z.x);
245                 y = fabs(z.y);
246                 if (x >= y)
247                 {
248                         r = y / x;
249                         w = sqrt(x) * sqrt(0.5 * (1.0 + sqrt(1.0 + r * r)));
250                 } else
251                 {
252                         r = x / y;
253                         w = sqrt(y) * sqrt(0.5 * (r + sqrt(1.0 + r * r)));
254                 }
255                 if (z.x >= 0.0)
256                 {
257                         c.x = w;
258                         c.y = z.y / (2.0 * w);
259                 } else
260                 {
261                         c.y = (z.y >= 0) ? w : -w;
262                         c.x = z.y / (2.0 * c.y);
263                 }
264                 return (c);
265         }
266 }
267 
268 Complex         cmplx_log(Complex z)
269 {
270         Complex         l = {Complex_id};
271 
272         l.x = 0.5 * log(z.x * z.x + z.y * z.y);
273         l.y = atan2(z.y, z.x);
274         return (l);
275 }
276 
277 Complex         cmplx_exp(Complex z)
278 {
279         Complex         l = {Complex_id};
280         double          expx;
281         double          tany;
282 
283         expx = exp(z.x);
284         tany = tan(z.y);
285         l.x = expx / sqrt(1.0 + tany * tany);
286         l.y = expx / sqrt(1.0 / (tany * tany) + 1.0);
287 
288         return (l);
289 }
290 
291 Complex         cmplx_sin(Complex z)
292 {
293         Complex         l = {Complex_id};
294 
295         l.x = sin(z.x) * cosh(z.y);
296         l.y = cos(z.x) * sinh(z.y);
297         return (l);
298 }
299 
300 Complex         cmplx_asin(Complex z)
301 {
302         Complex         l = {Complex_id};
303         double          alpha, beta, a, b;
304 
305         a = 0.5 * sqrt((z.x + 1.0) * (z.x + 1.0) + z.y * z.y);
306         b = 0.5 * sqrt((z.x - 1.0) * (z.x - 1.0) + z.y * z.y);
307         alpha = a + b;
308         beta = a - b;
309 
310         l.x = asin(beta);
311         l.y = log(alpha + sqrt(alpha * alpha - 1.0));
312         return (l);
313 }
314 
315 Complex         cmplx_cis(double theta)
316 {
317         return (cmplx(cos(theta), sin(theta)));
318 }
319 

~ [ source navigation ] ~ [ diff markup ] ~ [ identifier search ] ~ [ freetext search ] ~ [ file search ] ~

This page was automatically generated by the LXR engine.
Visit the LXR main site for more information.