From 40e5a5811c6ffce9b0974e93cdd927cbcf60c157 Mon Sep 17 00:00:00 2001 From: Joe Hunkeler Date: Tue, 11 Aug 2015 16:51:37 -0400 Subject: Repatch (from linux) of OSX IRAF --- unix/f2c/libf2c/c_div.c | 53 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 53 insertions(+) create mode 100644 unix/f2c/libf2c/c_div.c (limited to 'unix/f2c/libf2c/c_div.c') diff --git a/unix/f2c/libf2c/c_div.c b/unix/f2c/libf2c/c_div.c new file mode 100644 index 00000000..9463a43d --- /dev/null +++ b/unix/f2c/libf2c/c_div.c @@ -0,0 +1,53 @@ +#include "f2c.h" +#ifdef __cplusplus +extern "C" { +#endif + +#ifdef KR_headers +extern VOID sig_die(); +VOID c_div(c, a, b) +complex *a, *b, *c; +#else +extern void sig_die(const char*,int); +void c_div(complex *c, complex *a, complex *b) +#endif +{ + double ratio, den; + double abr, abi, cr; + + if( (abr = b->r) < 0.) + abr = - abr; + if( (abi = b->i) < 0.) + abi = - abi; + if( abr <= abi ) + { + if(abi == 0) { +#ifdef IEEE_COMPLEX_DIVIDE + float af, bf; + af = bf = abr; + if (a->i != 0 || a->r != 0) + af = 1.; + c->i = c->r = af / bf; + return; +#else + sig_die("complex division by zero", 1); +#endif + } + ratio = (double)b->r / b->i ; + den = b->i * (1 + ratio*ratio); + cr = (a->r*ratio + a->i) / den; + c->i = (a->i*ratio - a->r) / den; + } + + else + { + ratio = (double)b->i / b->r ; + den = b->r * (1 + ratio*ratio); + cr = (a->r + a->i*ratio) / den; + c->i = (a->i - a->r*ratio) / den; + } + c->r = cr; + } +#ifdef __cplusplus +} +#endif -- cgit