/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:32:04 */
/*FOR_C Options SET: ftn=u io=c no=p op=aimnv s=dbov str=l x=f - prototypes */
#include <math.h>
#include "fcrt.h"
#include "sfrenl.h"
#include <float.h>
#include <stdlib.h>
/*  File: sfrenl.[for|f|c]
 *  Contains procedures: SFRENC(), SFRENF(), SFRENG(), SFRENS()
 *  and private low-level procedure: SFREN1().
 * Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
 * ALL RIGHTS RESERVED.
 * Based on Government Sponsored Research NAS7-03001.
 *>> 2000-12-01 SFRENL Krogh Removed parameter PID2.
 *>> 2000-09-07 SFRENL Krogh Changed text on x to one on ax.
 *>> 1996-01-08 SFRENL WV Snyder Use SCSPXX for cos(Pi/2 x**2), etc.
 *>> 1995-11-03 SFRENL Krogh  Removed blanks in numbers for C conversion.
 *>> 1994-11-02 SFRENL Krogh  Changes to use M77CON
 *>> 1994-10-18 SFRENL WV Snyder More specializing instructions
 *>> 1993-02-25 SFRENL CLL. Edited to eliminate ENTRY and EQUIVALENCE.
 *>> 1992-09-15 SFRENL WV Snyder Specializing instructions
 *>> 1992-04-13 SFRENL WV Snyder Declare SFRENF, SFRENG, SFRENS
 *>> 1992-03-18 SFRENL WV Snyder Move declarations for coefficient arrays
 *>> 1992-01-24 SFRENL WV Snyder Original code
 *--S replaces "?": ?FRENC, ?FREN1, ?FRENF, ?FRENG, ?FRENS, ?FRENL
 *--&               ?CSPXX, ?SNPXX
 * Subprograms in this file compute the Fresnel Cosine and Sine
 * integrals C(x) and S(x), and the auxiliary functions f(x) and g(x),
 * for any X:
 *     SFRENC(X) for Fresnel integral C(X)
 *     SFRENS(X) for Fresnel integral S(x)
 *     SFRENF(X) for Fresnel integral auxiliary function f(x)
 *     SFRENG(X) for Fresnel integral auxiliary function g(x).
 *
 * Developed by W. V. Snyder, Jet Propulsion Laboratory, 24 January 1992.
 *
 * Ref: W. J. Cody, "Chebyshev Approximations for the Fresnel Integrals",
 * Mathematics of Computation, 1968, pp 450-453 plus Microfiche Suppl.
 * W. V. Snyder, "Algorithm 723: Fresnel Integrals," ACM Trans. Math.
 * Softw. 19, 4 (December 1993) 452-456.
 * Accuracies of highest order formulae, where E is relative error:
 *
 * Range           Function   -log10(E)   Function   -log10(E)
 * |X|<=1.2          C(x)       16.24       S(x)       17.26
 * 1.2<|X|<=1.6      C(x)       17.47       S(x)       18.66
 * 1.6<|X|<=1.9      f(x)       17.13       g(x)       16.25
 * 1.9<|X|<=2.4      f(x)       16.64       g(x)       15.65
 * 2.4<|X|           f(x)       16.89       g(x)       15.58
 *
 * Refer to Cody for accuracy of other approximations.
 *
 *     ================================================================== */
float /*FUNCTION*/ sfrenc(
float x)
{
	float sfrenc_v;
 
	sfrenc_v = sfren1( 1, x );
	return( sfrenc_v );
} /* end of function */
/*     ================================================================== */
float /*FUNCTION*/ sfrenf(
float x)
{
	float sfrenf_v;
 
	sfrenf_v = sfren1( 3, x );
	return( sfrenf_v );
} /* end of function */
/*     ================================================================== */
float /*FUNCTION*/ sfreng(
float x)
{
	float sfreng_v;
 
	sfreng_v = sfren1( 4, x );
	return( sfreng_v );
} /* end of function */
/*     ================================================================== */
float /*FUNCTION*/ sfrens(
float x)
{
	float sfrens_v;
 
	sfrens_v = sfren1( 2, x );
	return( sfrens_v );
} /* end of function */
/*     ================================================================== */
		/* PARAMETER translations */
#define	RPI	0.3183098861837906715377675267450287240689e0
#define	RPISQ	(RPI*RPI)
		/* end of PARAMETER translations */
 
float /*FUNCTION*/ sfren1(
long mode,
float x)
{
	LOGICAL32 wantc, wantf, wantg;
	long int _l0;
	float ax, cx, sfren1_v, sx, x4;
	static float largef, largeg, largex;
	static float bigx = -1.0e0;
	static float c = 0.0e0;
	static float f = 0.5e0;
	static float g = 0.5e0;
	static float s = 0.0e0;
	static float xsave = 0.0e0;
	static LOGICAL32 havec = TRUE;
	static LOGICAL32 havef = TRUE;
	static LOGICAL32 haveg = TRUE;
	static LOGICAL32 haves = TRUE;
	static LOGICAL32 needc[8]={TRUE,FALSE,TRUE,TRUE,TRUE,FALSE,FALSE,
	 FALSE};
	static LOGICAL32 needs[8]={FALSE,TRUE,TRUE,TRUE,FALSE,TRUE,FALSE,
	 FALSE};
	static LOGICAL32 needf[8]={FALSE,FALSE,TRUE,FALSE,TRUE,TRUE,TRUE,
	 FALSE};
	static LOGICAL32 needg[8]={FALSE,FALSE,FALSE,TRUE,TRUE,TRUE,FALSE,
	 TRUE};
	static float pc1[4-(0)+1]={9.999999999999999421e-1,-1.994608988261842706e-1,
	 1.761939525434914045e-2,-5.280796513726226960e-4,5.477113856826871660e-6};
	static float qc1[4-(1)+1]={4.727921120104532689e-2,1.099572150256418851e-3,
	 1.552378852769941331e-5,1.189389014228757184e-7};
	static float pc2[5-(0)+1]={1.00000000000111043640e0,-2.07073360335323894245e-1,
	 1.91870279431746926505e-2,-6.71376034694922109230e-4,1.02365435056105864908e-5,
	 -5.68293310121870728343e-8};
	static float qc2[5-(1)+1]={3.96667496952323433510e-2,7.88905245052359907842e-4,
	 1.01344630866749406081e-5,8.77945377892369265356e-8,4.41701374065009620393e-10};
	static float ps1[4-(0)+1]={5.2359877559829887021e-1,-7.0748991514452302596e-2,
	 3.8778212346368287939e-3,-8.4555728435277680591e-5,6.7174846662514086196e-7};
	static float qs1[4-(1)+1]={4.1122315114238422205e-2,8.1709194215213447204e-4,
	 9.6269087593903403370e-6,5.9528122767840998345e-8};
	static float ps2[5-(0)+1]={5.23598775598344165913e-1,-7.37766914010191323867e-2,
	 4.30730526504366510217e-3,-1.09540023911434994566e-4,1.28531043742724820610e-6,
	 -5.76765815593088804567e-9};
	static float qs2[5-(1)+1]={3.53398342767472162540e-2,6.18224620195473216538e-4,
	 6.87086265718620117905e-6,5.03090581246612375866e-8,2.05539124458579596075e-10};
	static float pf1[5-(0)+1]={3.1830975293580985290e-1,1.2226000551672961219e1,
	 1.2924886131901657025e2,4.3886367156695547655e2,4.1466722177958961672e2,
	 5.6771463664185116454e1};
	static float qf1[5-(1)+1]={3.8713003365583442831e1,4.1674359830705629745e2,
	 1.4740030733966610568e3,1.5371675584895759916e3,2.9113088788847831515e2};
	static float pf2[5-(0)+1]={3.183098818220169217e-1,1.958839410219691002e1,
	 3.398371349269842400e2,1.930076407867157531e3,3.091451615744296552e3,
	 7.177032493651399590e2};
	static float qf2[5-(1)+1]={6.184271381728873709e1,1.085350675006501251e3,
	 6.337471558511437898e3,1.093342489888087888e4,3.361216991805511494e3};
	static float pf3[6-(0)+1]={-9.675460329952532343e-2,-2.431275407194161683e1,
	 -1.947621998306889176e3,-6.059852197160773639e4,-7.076806952837779823e5,
	 -2.417656749061154155e6,-7.834914590078317336e5};
	static float qf3[6-(1)+1]={2.548289012949732752e2,2.099761536857815105e4,
	 6.924122509827708985e5,9.178823229918143780e6,4.292733255630186679e7,
	 4.803294784260528342e7};
	static float pg1[5-(0)+1]={1.013206188102747985e-1,4.445338275505123778e0,
	 5.311228134809894481e1,1.991828186789025318e2,1.962320379716626191e2,
	 2.054214324985006303e1};
	static float qg1[5-(1)+1]={4.539250196736893605e1,5.835905757164290666e2,
	 2.544731331818221034e3,3.481121478565452837e3,1.013794833960028555e3};
	static float pg2[5-(0)+1]={1.01321161761804586e-1,7.11205001789782823e0,
	 1.40959617911315524e2,9.08311749529593938e2,1.59268006085353864e3,
	 3.13330163068755950e2};
	static float qg2[5-(1)+1]={7.17128596939302198e1,1.49051922797329229e3,
	 1.06729678030580897e4,2.41315567213369742e4,1.15149832376260604e4};
	static float pg3[6-(0)+1]={-1.53989733819769316e-1,-4.31710157823357568e1,
	 -3.87754141746378493e3,-1.35678867813756347e5,-1.77758950838029676e6,
	 -6.66907061668636416e6,-1.72590224654836845e6};
	static float qg3[6-(1)+1]={2.86733194975899483e2,2.69183180396242536e4,
	 1.02878693056687506e6,1.62095600500231646e7,9.38695862531635179e7,
	 1.40622441123580005e8};
		/* OFFSET Vectors w/subscript range: 1 to dimension */
	LOGICAL32 *const Needc = &needc[0] - 1;
	LOGICAL32 *const Needf = &needf[0] - 1;
	LOGICAL32 *const Needg = &needg[0] - 1;
	LOGICAL32 *const Needs = &needs[0] - 1;
		/* end of OFFSET VECTORS */
 
	/*     MODE = 1 means compute C.
	 *     MODE = 2 means compute S.
	 *     MODE = 3 means compute F.
	 *     MODE = 4 means compute G.
	 *     ------------------------------------------------------------------
	 *                        Internal variables.
	 *
	 * RPI is the reciprocal of PI.
	 * RPISQ is the reciprocal of PI squared.
	 * AX is abs(x).
	 * BIGX is 1/sqrt(round-off).  If X > BIGX then to the working
	 *         precision x**2 is an integer (which we assume to be a multiple
	 *         of four), so cos(pi/2 * x**2) = 1, and sin(pi/2 * x**2) = 0.
	 * C and S are values of C(x) and S(x), respectively.
	 * CX and SX are cos(pi/2 * ax**2) and sin(pi/2 * ax**2), respectively.
	 * F and G are used to compute f(x) and g(x) when X > 1.6.
	 * HAVEC, HAVEF, HAVEG, HAVES are logical variables that indicate
	 *         whether the values stored in C, F, G and S correspond to the
	 *         value stored in X.  HAVEF indicates we have both F and G when
	 *         XSAVE .le. 1.6, and HAVEC indicates we have both C and S when
	 *         XSAVE .gt. 1.6.
	 * LARGEF is 1/(pi * underflow).  If X > LARGEF then f ~ 0.
	 * LARGEG is cbrt(1/(pi**2 * underflow)).  If X > LARGEG then g ~ 0.
	 * LARGEX is 1/sqrt(sqrt(underflow)).  If X > LARGEX then f ~ 1/(pi * x)
	 *         and g ~ 1/(pi**2 * x**3).
	 * MODE indicates the function to be computed: 1 = C(x), 2 = S(x),
	 *         3 = f(x), 4 = g(x).
	 * NEEDC, NEEDF, NEEDG, NEEDS are arrays indexed by MODE (MODE+4 when
	 *         X .gt. 1.6) that indicate what functions are needed.
	 * WANTC indicates whether C and S must be computed from F and G.
	 * WANTF and WANTG indicate we computed F and G on the present call.
	 * XSAVE is the most recently provided value of X.
	 * X4 is either X ** 4 or (1.0/X) ** 4.
	 *     If you change the order of approximation, you must change the
	 *     declarations and DATA statements for the coefficient arrays,
	 *     and the executable statements that evaluate the approximations.
	 *     ------------------------------------------------------------------ */
	/*++   default digits=16
	 *++   code for digits <= 3 is inactive
	 *      real             PC1(0:1), QC1(1:1)
	 *++   code for digits > 3 & digits <= 6 is inactive
	 *      real             PC1(0:2), QC1(1:2)
	 *++   code for digits > 6 & digits <= 11 is inactive
	 *      real             PC1(0:3), QC1(1:3)
	 *++   code for digits > 11 is active */
	/*++   end
	 *++   code for digits <= 2 is inactive
	 *      real             PC2(0:1), QC2(1:1)
	 *++   code for digits > 2 & digits <= 5 is inactive
	 *      real             PC2(0:2), QC2(1:2)
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      real             PC2(0:3), QC2(1:3)
	 *++   code for digits > 8 & digits <= 12 is inactive
	 *      real             PC2(0:4), QC2(1:4)
	 *++   code for digits > 12 is active */
	/*++   end
	 *++   code for digits <= 3 is inactive
	 *      real             PS1(0:1), QS1(1:1)
	 *++   code for digits > 3 & digits <= 7 is inactive
	 *      real             PS1(0:2), QS1(1:2)
	 *++   code for digits > 7 & digits <= 12 is inactive
	 *      real             PS1(0:3), QS1(1:3)
	 *++   code for digits > 12 is active */
	/*++   end
	 *++   code for digits <= 2 is inactive
	 *      real             PS2(0:1), QS2(1:1)
	 *++   code for digits > 2 & digits <= 5 is inactive
	 *      real             PS2(0:2), QS2(1:2)
	 *++   code for digits > 5 & digits <= 9 is inactive
	 *      real             PS2(0:3), QS2(1:3)
	 *++   code for digits > 9 & digits <= 14 is inactive
	 *      real             PS2(0:4), QS2(1:4)
	 *++   code for digits > 14 is active */
	/*++   end
	 *++   code for digits <= 5 is inactive
	 *      real             PF1(0:1), QF1(1:1)
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      real             PF1(0:2), QF1(1:2)
	 *++   code for digits > 8 & digits <= 11 is inactive
	 *      real             PF1(0:3), QF1(1:3)
	 *++   code for digits > 11 & digits <= 14 is inactive
	 *      real             PF1(0:4), QF1(1:4)
	 *++   code for digits > 14 is active */
	/*++   end
	 *++   code for digits <= 5 is inactive
	 *      real             PF2(0:1), QF2(1:1)
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      real             PF2(0:2), QF2(1:2)
	 *++   code for digits > 8 & digits <= 11 is inactive
	 *      real             PF2(0:3), QF2(1:3)
	 *++   code for digits > 11 & digits <= 14 is inactive
	 *      real             PF2(0:4), QF2(1:4)
	 *++   code for digits > 14 is active */
	/*++   end
	 *++   code for digits <= 3 is inactive
	 *      real             PF3(0:0)
	 *++   code for digits > 3 & digits <= 6 is inactive
	 *      real             PF3(0:1), QF3(1:1)
	 *++   code for digits > 6 & digits <= 9 is inactive
	 *      real             PF3(0:2), QF3(1:2)
	 *++   code for digits > 9 & digits <= 11 is inactive
	 *      real             PF3(0:3), QF3(1:3)
	 *++   code for digits > 11 & digits <= 13 is inactive
	 *      real             PF3(0:4), QF3(1:4)
	 *++   code for digits > 13 & digits <= 15 is inactive
	 *      real             PF3(0:5), QF3(1:5)
	 *++   code for digits > 15 is active */
	/*++   end
	 *++   code for digits <= 4 is inactive
	 *      real             PG1(0:1), QG1(1:1)
	 *++   code for digits > 4 & digits <= 7 is inactive
	 *      real             PG1(0:2), QG1(1:2)
	 *++   code for digits > 7 & digits <= 10 is inactive
	 *      real             PG1(0:3), QG1(1:3)
	 *++   code for digits > 10 & digits <= 13 is inactive
	 *     real             PG1(0:4), QG1(1:4)
	 *++   code for digits > 13 is active */
	/*++   end
	 *++   code for digits <= 4 is inactive
	 *      real             PG2(0:1), QG2(1:1)
	 *++   code for digits > 4 & digits <= 7 is inactive
	 *      real             PG2(0:2), QG2(1:2)
	 *++   code for digits > 7 & digits <= 10 is inactive
	 *      real             PG2(0:3), QG2(1:3)
	 *++   code for digits > 10 & digits <= 13 is inactive
	 *      real             PG2(0:4), QG2(1:4)
	 *++   code for digits > 13 is active */
	/*++   end
	 *++   code for digits <= 3 is inactive
	 *      real             PG3(0:0)
	 *++   code for digits > 3 & digits <= 5 is inactive
	 *      real             PG3(0:1), QG3(r:1)
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      real             PG3(0:2), QG3(1:2)
	 *++   code for digits > 8 & digits <= 10 is inactive
	 *      real             PG3(0:3), QG3(1:3)
	 *++   code for digits > 10 & digits <= 12 is inactive
	 *      real             PG3(0:4), QG3(1:4)
	 *++   code for digits = 13 is inactive
	 *      real             PG3(0:5), QG3(1:5)
	 *++   code for digits > 13 is active */
	/*++   end
	 * */
	/*        C(x)    S(x)    f(x)    g(x)    C(x)    S(x)    f(x)    g(x) */
 
	/*     Coefficients for C(x), |X| <= 1.2
	 *++   code for digits <= 3 is inactive
	 *      data pc1 / 1.00053e0, -1.12353e-1/, qc1 /1.38937e-1/
	 *++   code for digits > 3 & digits <= 6 is inactive
	 *      data pc1 / 9.99999896e-1, -1.63090954e-1, 1.06388604e-2/
	 *      data qc1 / 8.36467414e-2,  3.10155884e-3/
	 *++   code for digits > 6 & digits <= 11 is inactive
	 *      data pc1 / 1.0000000000042e0,  -1.8651127631106e-1,
	 *     *           1.5065663274457e-2, -3.1058074693185e-4/
	 *      data qc1 / 6.0228833908128e-2,  1.7410300558198e-3,
	 *     *            .6308617899210e-5/
	 *++   code for digits > 11 is active */
	/*++   end
	 *
	 *     Coefficients for C(x), 1.2 < |X| <= 1.6
	 *++   code for digits <= 2 is inactive
	 *      data pc2 / 1.3139e0, -5.8827e-2/, qc2 / 4.7324e-1/
	 *++   code for digits > 2 & digits <= 5 is inactive
	 *      data pc2 / 9.9790890e-1, -1.5391895e-1, 9.9817933e-3/
	 *      data qc2 / 8.9878647e-2,  5.6003071e-3/
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      data pc2 / 1.00000440109e0, -1.83413851908e-1,
	 *     *           1.45776170828e-2,-2.78980705270e-4/
	 *      data qc2 / 6.33347445637e-2, 2.01245738369e-3,
	 *     *           4.03949004646e-5/
	 *++   code for digits > 8 & digits <= 12 is inactive
	 *      data pc2 / 9.999999966496876e-1, -1.980300987022688e-1,
	 *     *           1.735518748450023e-2, -5.069442297935788e-4,
	 *     *           5.059962254678234e-6/
	 *      data qc2 / 4.871000308997918e-2,  1.188406974004084e-3,
	 *     *           1.824527635843850e-5,  1.678314257874103e-7/
	 *++   code for digits > 12 is active */
	/*++   end
	 *
	 *     Coefficients for S(x), |X| <= 1.2
	 *++   code for digits <= 3 is inactive
	 *      data ps1 / 5.23677e-1, -4.67900e-2/, qs1 / 8.81622e-2/
	 *++   code for digits > 3 & digits <= 7 is inactive
	 *      data ps1 / 5.235987665e-1,-5.837763961e-2, 2.165920196e-3/
	 *      data qs1 / 6.474944765e-2, 1.713287588e-3/
	 *++   code for digits > 7 & digits <= 12 is inactive
	 *      data ps1 / 5.23598775598566e-1, -6.59149581139046e-2,
	 *     *           3.21501649828293e-3, -4.88704436240178e-5/
	 *      data qs1 / 5.03546388670085e-2,  1.17835980356588e-3,
	 *     *           1.37089875826980e-5/
	 *++   code for digits > 12 is active */
	/*++   end
	 *
	 *     coefficients for S(x), 1.2 < |X| <= 1.6
	 *++   code for digits <= 2 is inactive
	 *      data ps2 / 5.4766e-1, -3.7151e-2/, qs2 / 1.4559e-1/
	 *++   code for digits > 2 & digits <= 5 is inactive
	 *      data ps2 / 5.2343994e-1, -5.4828347e-2, 1.9020881e-3/
	 *      data qs2 / 7.1104704e-2,  2.5596274e-3/
	 *++   code for digits > 5 & digits <= 9 is inactive
	 *      data ps2 / 5.23599040498e-1, -6.46593392426e-2,
	 *     *           3.08030794361e-3, -4.40800854418e-5/
	 *      data qs2 / 5.27536681685e-2,  1.34311026821e-3,
	 *     *           1.90476612849e-5/
	 *++   code for digits > 9 & digits <= 14 is inactive
	 *      data ps2 / 5.2359877543509178e-1,-7.0149076634833662e-2,
	 *     *           3.8031581605987038e-3,-8.0964948714408156e-5,
	 *     *           6.1908080210052772e-7/
	 *      data qs2 / 4.2268067370395487e-2, 8.7642753831073237e-4,
	 *     *           1.1088542889789282e-5, 7.8725829545478464e-8/
	 *++   code for digits > 14 is active */
	/*++   end
	 *
	 *     coefficients for f(x), 1.6 < |X| <= 1.9
	 *++   code for digits <= 5 is inactive
	 *      data pf1 / 3.1803519e-1, 4.2352332e-1/, qf1 / 1.6015403e0/
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      data pf1 / 3.18285252094e-1,  2.02860275713e0,
	 *     *           1.08108136048e0 /
	 *      data qf1 / 6.67171548135e0,   4.52997553972e0 /
	 *++   code for digits > 8 & digits <= 11 is inactive
	 *      data pf1 / 3.1830646311448e-1,  4.6077249266684e0,
	 *     *           1.1914578438074e1,   3.5555073665830e0 /
	 *      data qf1 / 1.4778464938936e1,   4.0902971163705e1,
	 *     *           1.6130432784819e1 /
	 *++   code for digits > 11 & digits <= 14 is inactive
	 *      data pf1 / 3.1830926850490599e-1,  8.0358812280394156e0,
	 *     *           4.8034065557792487e1,   6.9853426160102065e1,
	 *     *           1.3530423554038784e1 /
	 *      data qf1 / 2.5549161843579529e1,   1.5761100558012250e2,
	 *     *           2.4956199380517229e2,   6.5563064008391554e1 /
	 *++   code for digits > 14 is active */
	/*++   end
	 *
	 *     coefficients for f(x), 1.9 < |X| <= 2.4
	 *++   code for digits <= 5 is inactive
	 *      data pf2 / 3.1825112e-1, 5.8395951e-1/, qf2 / 2.1243944e0 /
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      data pf2 / 3.1830699932e-1, 2.9993457087e0, 2.3956340010e0/
	 *      data qf2 / 9.7254517756e0,  9.4814077696e0 /
	 *++   code for digits > 8 & digits <= 11 is inactive
	 *      data pf2 / 3.1830963944521e-1,  7.0770431878327e0,
	 *     *           2.8756083262903e1,   1.3583685742326e1 /
	 *      data qf2 / 2.2536994405207e1,   9.6127568469278e1,
	 *     *           5.7407028004031e1 /
	 *++   code for digits > 11 & digits <= 14 is inactive
	 *      data pf2 / 3.183098568640159e-1,  1.265129469683175e1,
	 *     *           1.219467616498339e2,   2.910033655512762e2,
	 *     *           9.278397828631516e1/
	 *      data qf2 / 4.004915302781009e1,   3.942059697951583e2,
	 *     *           1.001368403495691e3,   4.142676224222433e2/
	 *++   code for digits > 14 is active */
	/*++   end
	 *
	 *     coefficients for f(x), 2.4 < |X|
	 *++   code for digits <= 3 is inactive
	 *      data pf3 /-8.97969e-2/
	 *++   code for digits > 3 & digits <= 6 is inactive
	 *      data pf3 / -9.67122165e-2, -3.73565920e-1/
	 *      data qf3 / 7.29466572e0 /
	 *++   code for digits > 6 & digits <= 9 is inactive
	 *      data pf3 /-9.67541443648e-2, -2.26566293818e0,
	 *      *         -3.53429821084e0 /
	 *      data qf3 / 2.69596939977e1,   9.72883957469e1 /
	 *++   code for digits > 9 & digits <= 11 is inactive
	 *      data pf3 /-9.6754596389025e-2, -5.5254943840897e0,
	 *     *          -5.8606282086171e1,  -5.2235560918394e1 /
	 *      data qf3 / 6.0654483020979e1,   7.8528841711294e2,
	 *     *           1.8592498578831e3 /
	 *++   code for digits > 11 & digits <= 13 is inactive
	 *      data pf3 /-9.675460316952504e-2, -1.023876428129288e1,
	 *     *          -2.712579634037998e2,  -1.766119345282127e3,
	 *     *          -1.043464426656267e3 /
	 *      data qf3 / 1.093682244053459e2,   3.155843461920557e3,
	 *     *           2.625634316204420e4,   4.578252057246393e4 /
	 *++   code for digits > 13 & digits <= 15 is inactive
	 *      data pf3 /-9.67546032967090380e-2,-1.64797712841245767e1,
	 *     *          -8.16343401784374598e2, -1.34922028171857248e4,
	 *     *          -6.13547113614699772e4, -2.61294753225141779e4/
	 *      data qf3 / 1.73871690673649114e2,  9.01827596231524147e3,
	 *     *           1.65946462621853184e5,  1.00105478900791339e6,
	 *     *           1.37012364817225972e6 /
	 *++   code for digits > 15 is active */
	/*++   end
	 *
	 *     coefficients for g(x), 1.6 < |X| <= 1.9
	 *++   code for digits <= 4 is inactive
	 *      data pg1 / 1.007011e-1, 1.457780e-1/, qg1 / 2.700253e0 /
	 *++   code for digits > 4 & digits <= 7 is inactive
	 *      data pg1 / 1.012500483e-1, 7.735207446e-1, 3.878282770e-1/
	 *      data qg1 / 9.101956178e0,  1.002234185e1 /
	 *++   code for digits > 7 & digits <= 10 is inactive
	 *      data pg1 / 1.013095436817e-1,  1.731937984173e0,
	 *     *           5.036588245265e0,   1.290707246507e0 /
	 *      data qg1 / 1.860077430076e1,   6.901951935725e1,
	 *     *           4.308918659989e1 /
	 *++   code for digits > 10 & digits <= 13 is inactive
	 *      data pg1 / 1.013188096509180e-1,  2.966220554725899e0,
	 *     *           2.015435299505393e1,   3.155605679387908e1,
	 *     *           4.916012830646366e0 /
	 *      data qg1 / 3.079178736724045e1,   2.362874318980473e2,
	 *     *           4.953861258248338e2,   2.026144493403599e2 /
	 *++   code for digits > 13 is active */
	/*++   end
	 *
	 *     coefficients for g(x), 1.9 < |X| <= 2.4
	 *++   code for digits <= 4 is inactive
	 *      data pg2 / 1.011711e-1,  2.239630e-1/, qg2 / 3.609237e0 /
	 *++   code for digits > 4 & digits <= 7 is inactive
	 *      data pg2 / 1.013115463e-1,  1.182021191e0,  9.886315969e-1/
	 *      data qg2 / 1.317219285e1,   2.101104994e1 /
	 *++   code for digits > 7 & digits <= 10 is inactive
	 *      data pg2 / 1.013202025653e-1,  2.690767013770e0,
	 *     *           1.283624749271e1,   5.791337587723e0 /
	 *      data qg2 / 2.807457005500e1,   1.598704313522e2,
	 *     *           1.525630385045e2 /
	 *++   code for digits > 10 & digits <= 13 is inactive
	 *      data pg2 / 1.013210509409046e-1,  4.682769769757399e0,
	 *     *           5.241351613346472e1,   1.411735944550041e2,
	 *     *           4.019756012788710e1 /
	 *      data qg2 / 4.773652966704634e1,   5.802036967947208e2,
	 *     *           1.947179327244406e3,   1.266041236960445e3 /
	 *++   code for digits > 13 is active */
	/*++   end
	 *
	 *     coefficients for g(x), 2.4 < |X|
	 *++   code for digits <= 3 is inactive
	 *      data pg3 /-1.3549e-1/
	 *++   code for digits > 3 & digits <= 5 is inactive
	 *      data pg3 /-1.5382824e-1, -6.4547464e-1/, qg3 / 1.0290328e1/
	 *++   code for digits > 5 & digits <= 8 is inactive
	 *      data pg3 /-1.5398760302e-1, -4.2772857997e0,
	 *     *          -6.5940181141e0 /
	 *      data qg3 / 3.4149986477e1,   1.7071708816e2 /
	 *++   code for digits > 8 & digits <= 10 is inactive
	 *      data pg3 /-1.539896971616e-1, -1.024638314446e1,
	 *     *          -1.252507402120e2,  -1.029189676144e2 /
	 *      data qg3 / 7.292229303754e1,   1.186528673248e3,
	 *     *           3.844430908476e3 /
	 *++   code for digits > 10 & digits <= 12 is inactive
	 *      data pg3 /-1.53989733057971e-1, -1.86483223831639e1,
	 *     *          -5.66882778026550e2,  -4.16714647017489e3,
	 *     *          -2.14678074364341e3 /
	 *      data qg3 / 1.27484298075273e2,   4.40258858159927e3,
	 *     *           4.57583830694684e4,   1.08207883281275e5 /
	 *++   code for digits = 13 is inactive
	 *      data pg3 /-1.539897338019247e-1, -2.959078231855258e1,
	 *     *          -1.661867087183632e3,  -3.112796454920657e4,
	 *     *          -1.573536286819766e5,  -5.574884113746041e4 /
	 *      data qg3 / 1.985439813544440e2,   1.196693134508007e4,
	 *     *           2.625573864512749e5,   1.969055829311119e6,
	 *     *           3.629081333313312e6 /
	 *++   code for digits > 13 is active */
	/*++   end
	 *     ------------------------------------------------------------------ */
	if (bigx < 0.0e0)
	{
		bigx = 1.0e0/FLT_EPSILON;
		largef = RPI/FLT_MIN;
		largeg = powf(RPI*largef,1.0e0/3.0e0);
		largex = 1.0e0/sqrtf( sqrtf( FLT_MIN ) );
	}
	if (x != xsave)
	{
		havec = FALSE;
		havef = FALSE;
		haveg = FALSE;
		haves = FALSE;
	}
	ax = fabsf( x );
	if (ax <= 1.6e0)
	{
		x4 = powif(ax,4);
		if (Needc[mode] && !havec)
		{
			if (ax <= 1.2e0)
			{
				/*++             code for digits <= 3 is inactive
				 *               c = x * (pc1(1)*x4+pc1(0)) / (qc1(1)*x4+1.0e0)
				 *++             code for digits > 3 & digits <= 6 is inactive
				 *               c = x * ((pc1(2)*x4+pc1(1))*x4+pc1(0))
				 *     1               / ((qc1(2)*x4+qc1(1))*x4+1.0e0)
				 *++             code for digits > 6 & digits <= 11 is inactive
				 *               c = x * (((pc1(3)*x4+pc1(2))*x4+pc1(1))*x4+pc1(0))
				 *     1               / (((qc1(3)*x4+qc1(2))*x4+qc1(1))*x4+1.0e0)
				 *++             code for digits > 11 is active */
				c = x*((((pc1[4]*x4 + pc1[3])*x4 + pc1[2])*x4 + pc1[1])*
				 x4 + pc1[0])/((((qc1[3]*x4 + qc1[2])*x4 + qc1[1])*
				 x4 + qc1[0])*x4 + 1.0e0);
				/*++             end */
			}
			else
			{
				/*++             code for digits <= 2 is inactive
				 *               c = x * (pc2(1)*x4+pc2(0)) / (qc2(1)*x4+1.0e0)
				 *++             code for digits > 2 & digits <= 5 is inactive
				 *               c = x * ((pc2(2)*x4+pc2(1))*x4+pc2(0))
				 *     1               / ((qc2(2)*x4+qc2(1))*x4+1.0e0)
				 *++             code for digits > 5 & digits <= 8 is inactive
				 *               c = x * (((pc2(3)*x4+pc2(2))*x4+pc2(1))*x4+pc2(0))
				 *     1               / (((qc2(3)*x4+qc2(2))*x4+qc2(1))*x4+1.0e0)
				 *++             code for digits > 8 & digits <= 12 is inactive
				 *               c = x * ((((pc2(4)*x4+pc2(3))*x4+pc2(2))*x4+pc2(1))*x4+
				 *     1                   pc2(0))
				 *     2               / ((((qc2(4)*x4+qc2(3))*x4+qc2(2))*x4+qc2(1))*x4+
				 *     3                   1.0e0)
				 *++             code for digits > 12 is active */
				c = x*(((((pc2[5]*x4 + pc2[4])*x4 + pc2[3])*x4 + pc2[2])*
				 x4 + pc2[1])*x4 + pc2[0])/(((((qc2[4]*x4 + qc2[3])*
				 x4 + qc2[2])*x4 + qc2[1])*x4 + qc2[0])*x4 + 1.0e0);
				/*++             end */
			}
			havec = TRUE;
		}
		if (Needs[mode] && !haves)
		{
			if (ax <= 1.2e0)
			{
				/*++             code for digits <= 3 is inactive
				 *               s = x**3*(ps1(1)*x4+ps1(0)) / (qs1(1)*x4+1.0e0)
				 *++             code for digits > 3 & digits <= 7 is inactive
				 *               s = x**3*((ps1(2)*x4+ps1(1))*x4+ps1(0))
				 *     1                / ((qs1(2)*x4+qs1(1))*x4+1.0e0)
				 *++             code for digits > 7 & digits <= 12 is inactive
				 *               s = x**3*(((ps1(3)*x4+ps1(2))*x4+ps1(1))*x4+ps1(0))
				 *     1                / (((qs1(3)*x4+qs1(2))*x4+qs1(1))*x4+1.0e0)
				 *++             code for digits > 12 is active */
				s = CUBE(x)*((((ps1[4]*x4 + ps1[3])*x4 + ps1[2])*x4 +
				 ps1[1])*x4 + ps1[0])/((((qs1[3]*x4 + qs1[2])*x4 +
				 qs1[1])*x4 + qs1[0])*x4 + 1.0e0);
				/*++             end */
			}
			else
			{
				/*++             code for digits <= 2 is inactive
				 *               s = x**3*(ps2(1)*x4+ps2(0)) / (qs2(1)*x4+1.0e0)
				 *++             code for digits > 2 & digits <= 5 is inactive
				 *               s = x**3*((ps2(2)*x4+ps2(1))*x4+ps2(0))
				 *     1              /   ((qs2(2)*x4+qs2(1))*x4+1.0e0)
				 *++             code for digits > 5 & digits <= 9 is inactive
				 *               s = x**3*(((ps2(3)*x4+ps2(2))*x4+ps2(1))*x4+ps2(0))
				 *     1              /   (((qs2(3)*x4+qs2(2))*x4+qs2(1))*x4+1.0e0)
				 *++             code for digits > 9 & digits <= 14 is inactive
				 *               s = x**3*((((ps2(4)*x4+ps2(3))*x4+ps2(2))*x4+ps2(1))*x4+
				 *     1                      ps2(0))
				 *     2              /   ((((qs2(4)*x4+qs2(3))*x4+qs2(2))*x4+qs2(1))*x4+
				 *     3                      1.0e0)
				 *++             code for digits > 14 is active */
				s = CUBE(x)*(((((ps2[5]*x4 + ps2[4])*x4 + ps2[3])*
				 x4 + ps2[2])*x4 + ps2[1])*x4 + ps2[0])/(((((qs2[4]*
				 x4 + qs2[3])*x4 + qs2[2])*x4 + qs2[1])*x4 + qs2[0])*
				 x4 + 1.0e0);
				/*++             end */
			}
			haves = TRUE;
		}
		if ((Needf[mode] || Needg[mode]) && !havef)
		{
			cx = scspxx( ax );
			sx = ssnpxx( ax );
			f = (0.5e0 - s)*cx - (0.5e0 - c)*sx;
			g = (0.5e0 - c)*cx + (0.5e0 - s)*sx;
			havef = TRUE;
		}
	}
	else
	{
		if (ax <= largex)
		{
			x4 = powif(1.0e0/ax,4);
			wantf = Needf[mode + 4] && !havef;
			if (wantf)
			{
				if (ax <= 1.9e0)
				{
					/*++                code for digits <= 5 is inactive
					 *                  f = (pf1(1)*x4+pf1(0)) / ((qf1(1)*x4+1.0e0) * ax)
					 *++                code for digits > 5 & digits <= 8 is inactive
					 *                  f = ((pf1(2)*x4+pf1(1))*x4+pf1(0))
					 *     1             / (((qf1(2)*x4+qf1(1))*x4+1.0e0) * ax)
					 *++                code for digits > 8 & digits <= 11 is inactive
					 *                  f = (((pf1(3)*x4+pf1(2))*x4+pf1(1))*x4+pf1(0))
					 *     1             / ((((qf1(3)*x4+qf1(2))*x4+qf1(1))*x4+1.0e0) * ax)
					 *++                code for digits > 11 & digits <= 14 is inactive
					 *                  f = ((((pf1(4)*x4+pf1(3))*x4+pf1(2))*x4+pf1(1))*x4+
					 *     1                    pf1(0))
					 *     2             / (((((qf1(4)*x4+qf1(3))*x4+qf1(2))*x4+qf1(1))*x4+
					 *     3                    1.0e0) * ax)
					 *++                code for digits > 14 is active */
					f = (((((pf1[5]*x4 + pf1[4])*x4 + pf1[3])*x4 +
					 pf1[2])*x4 + pf1[1])*x4 + pf1[0])/((((((qf1[4]*
					 x4 + qf1[3])*x4 + qf1[2])*x4 + qf1[1])*x4 + qf1[0])*
					 x4 + 1.0e0)*ax);
					/*++                end */
				}
				else if (ax <= 2.4)
				{
					/*++                code for digits <= 5 is inactive
					 *                  f = (pf2(1)*x4+pf2(0)) / ((qf2(1))*x4+1.0e0) * ax)
					 *++                code for digits > 5 & digits <= 8 is inactive
					 *                  f = ((pf2(2)*x4+pf2(1))*x4+pf2(0))
					 *     1             / (((qf2(2)*x4+qf2(1))*x4+1.0e0) * ax)
					 *++                code for digits > 8 & digits <= 11 is inactive
					 *                  f = (((pf2(3)*x4+pf2(2))*x4+pf2(1))*x4+pf2(0))
					 *     1             / ((((qf2(3)*x4+qf2(2))*x4+qf2(1))*x4+1.0e0) * ax)
					 *++                code for digits > 11 & digits <= 14 is inactive
					 *                  f = ((((pf2(4)*x4+pf2(3))*x4+pf2(2))*x4+pf2(1))*x4+
					 *     1                    pf2(0))
					 *     2             / (((((qf2(4)*x4+qf2(3))*x4+qf2(2))*x4+qf2(1))*x4+
					 *     3                    1.0e0) * ax)
					 *++                code for digits > 14 is active */
					f = (((((pf2[5]*x4 + pf2[4])*x4 + pf2[3])*x4 +
					 pf2[2])*x4 + pf2[1])*x4 + pf2[0])/((((((qf2[4]*
					 x4 + qf2[3])*x4 + qf2[2])*x4 + qf2[1])*x4 + qf2[0])*
					 x4 + 1.0e0)*ax);
					/*++                end */
				}
				else
				{
					/*++                code for digits <= 3 is inactive
					 *                  f = (rpi + x4*pf3(0)) / ax
					 *++                code for digits > 3 & digits <= 6 is inactive
					 *                  f = (rpi + x4*(pf3(1)*x4+pf3(0)) / (qf3(1)*x4+1.0e0))
					 *                      / ax
					 *++                code for digits > 6 & digits <= 9 is inactive
					 *                  f = (rpi +
					 *     1              x4*((pf3(2)*x4+pf3(1))*x4+pf3(0))
					 *     2            /    ((qf3(2)*x4+qf3(1))*x4+1.0e0)) / ax
					 *++                code for digits > 9 & digits <= 11 is inactive
					 *                  f = (rpi +
					 *     1              x4*(((pf3(3)*x4+pf3(2))*x4+pf3(1))*x4+pf3(0))
					 *     2            /    (((qf3(3)*x4+qf3(2))*x4+qf3(1))*x4+1.0e0)) / ax
					 *++                code for digits > 11 & digits <= 13 is inactive
					 *                  f = (rpi +
					 *     1              x4*((((pf3(4)*x4+pf3(3))*x4+pf3(2))*x4+pf3(1))*x4+
					 *     2                   pf3(0))
					 *     3            /    ((((qf3(4)*x4+qf3(3))*x4+qf3(2))*x4+qf3(1))*x4+
					 *     4                   1.0e0)) / ax
					 *++                code for digits > 13 & digits <= 15 is inactive
					 *                  f = (rpi +
					 *     1              x4*(((((pf3(5)*x4+pf3(4))*x4+pf3(3))*x4+pf3(2))*x4+
					 *     2                   pf3(1))*x4+pf3(0))
					 *     3            /    (((((qf3(5)*x4+qf3(4))*x4+qf3(3))*x4+qf3(2))*x4+
					 *     4                   qf3(1))*x4+1.0e0)) / ax
					 *++                code for digits > 15 is active */
					f = (RPI + x4*((((((pf3[6]*x4 + pf3[5])*x4 + pf3[4])*
					 x4 + pf3[3])*x4 + pf3[2])*x4 + pf3[1])*x4 + pf3[0])/
					 ((((((qf3[5]*x4 + qf3[4])*x4 + qf3[3])*x4 + qf3[2])*
					 x4 + qf3[1])*x4 + qf3[0])*x4 + 1.0e0))/ax;
					/*++                end */
				}
				havef = TRUE;
			}
			wantg = Needg[mode + 4] && !haveg;
			if (wantg)
			{
				if (ax <= 1.9e0)
				{
					/*++                code for digits <=4 is inactive
					 *                  g = (pg1(1)*x4+pg1(0)) / ((qg1(1)*x4+1.0e0) * ax**3)
					 *++                code for digits > 4 & digits <= 7 is inactive
					 *                  g = ((pg1(2)*x4+pg1(1))*x4+pg1(0))
					 *     1             / (((qg1(2)*x4+qg1(1))*x4+1.0e0) * ax**3)
					 *++                code for digits > 7 & digits <= 10 is inactive
					 *                  g = (((pg1(3)*x4+pg1(2))*x4+pg1(1))*x4+pg1(0))
					 *     1             / ((((qg1(3)*x4+qg1(2))*x4+qg1(1))*x4+1.0e0) *ax**3)
					 *++                code for digits > 10 & digits <= 13 is inactive
					 *                  g = ((((pg1(4)*x4+pg1(3))*x4+pg1(2))*x4+pg1(1))*x4+
					 *     1                    pg1(0))
					 *     2             / (((((qg1(4)*x4+qg1(3))*x4+qg1(2))*x4+qg1(1))*x4+
					 *     3                    1.0e0) * ax**3)
					 *++                code for digits > 13 is active */
					g = (((((pg1[5]*x4 + pg1[4])*x4 + pg1[3])*x4 +
					 pg1[2])*x4 + pg1[1])*x4 + pg1[0])/((((((qg1[4]*
					 x4 + qg1[3])*x4 + qg1[2])*x4 + qg1[1])*x4 + qg1[0])*
					 x4 + 1.0e0)*CUBE(ax));
					/*++                end */
				}
				else if (ax <= 2.4e0)
				{
					/*++                code for digits <=4 is inactive
					 *                  g = (pg2(1)*x4+pg2(0)) / ((qg2(1)*x4+1.0e0) * ax**3)
					 *++                code for digits > 4 & digits <= 7 is inactive
					 *                  g = ((pg2(2)*x4+pg2(1))*x4+pg2(0))
					 *     1             / (((qg2(2)*x4+qg2(1))*x4+1.0e0) * ax**3)
					 *++                code for digits > 7 & digits <= 10 is inactive
					 *                  g = (((pg2(3)*x4+pg2(2))*x4+pg2(1))*x4+pg2(0))
					 *     1             / ((((qg2(3)*x4+qg2(2))*x4+qg2(1))*x4+1.0e0) *ax**3)
					 *++                code for digits > 10 & digits <= 13 is inactive
					 *                  g = ((((pg2(4)*x4+pg2(3))*x4+pg2(2))*x4+pg2(1))*x4+
					 *     1                    pg2(0))
					 *     2             / (((((qg2(4)*x4+qg2(3))*x4+qg2(2))*x4+qg2(1))*x4+
					 *     3                    1.0e0) * ax**3)
					 *++                code for digits > 13 is active */
					g = (((((pg2[5]*x4 + pg2[4])*x4 + pg2[3])*x4 +
					 pg2[2])*x4 + pg2[1])*x4 + pg2[0])/((((((qg2[4]*
					 x4 + qg2[3])*x4 + qg2[2])*x4 + qg2[1])*x4 + qg2[0])*
					 x4 + 1.0e0)*CUBE(ax));
					/*++                end */
				}
				else
				{
					/*++                code for digits <=3 is inactive
					 *                  g = (rpisq + x4*pg3(0)) / ax**3
					 *++                code for digits > 3 & digits <= 5 is inactive
					 *                  g = (rpisq + x4*(pg3(1)*x4+pg3(0))
					 *     1            /    (qg3(1)*x4+1.0e0)) / ax**3
					 *++                code for digits > 5 & digits <= 8 is inactive
					 *                  g = (rpisq + x4*((pg3(2)*x4+pg3(1))*x4+pg3(0))
					 *     1            /    ((qg3(2)*x4+qg3(1))*x4+1.0e0)) / ax**3
					 *++                code for digits > 8 & digits <= 10 is inactive
					 *                  g = (rpisq +
					 *     1              x4*(((pg3(3)*x4+pg3(2))*x4+pg3(1))*x4+pg3(0))
					 *     2            /    (((qg3(3)*x4+qg3(2))*x4+qg3(1))*x4+1.0e0))
					 *     3            / ax**3
					 *++                code for digits > 10 & digits <= 12 is inactive
					 *                  g = (rpisq +
					 *     1              x4*((((pg3(4)*x4+pg3(3))*x4+pg3(2))*x4+pg3(1))*x4+
					 *     2                     pg3(0))
					 *     3            /    ((((qg3(4)*x4+qg3(3))*x4+qg3(2))*x4+qg3(1))*x4+
					 *     4                     1.0e0)) / ax**3
					 *++                code for digits = 13 is inactive
					 *                  g = (rpisq +
					 *     1              x4*(((((pg3(5)*x4+pg3(4))*x4+pg3(3))*x4+pg3(2))*x4+
					 *     2                      pg3(1))*x4+pg3(0))
					 *     3            /    (((((qg3(5)*x4+qg3(4))*x4+qg3(3))*x4+qg3(2))*x4+
					 *     4                      qg3(1))*x4+1.0e0)) / ax**3
					 *++                code for digits > 13 is active */
					g = (RPISQ + x4*((((((pg3[6]*x4 + pg3[5])*x4 +
					 pg3[4])*x4 + pg3[3])*x4 + pg3[2])*x4 + pg3[1])*
					 x4 + pg3[0])/((((((qg3[5]*x4 + qg3[4])*x4 + qg3[3])*
					 x4 + qg3[2])*x4 + qg3[1])*x4 + qg3[0])*x4 + 1.0e0))/
					 CUBE(ax);
					/*++                end */
				}
				haveg = TRUE;
			}
		}
		else
		{
			wantf = Needf[mode];
			if (wantf)
			{
				if (x <= largef)
				{
					f = RPI/ax;
				}
				else
				{
					f = 0.0e0;
				}
			}
			wantg = Needg[mode];
			if (wantg)
			{
				if (x <= largeg)
				{
					g = RPISQ/CUBE(ax);
				}
				else
				{
					g = 0.0e0;
				}
			}
		}
		wantc = (Needc[mode + 4] || Needs[mode + 4]) && !havec;
		if (wantc || x < 0.0e0)
		{
			cx = scspxx( ax );
			sx = ssnpxx( ax );
			if (wantc)
			{
				c = 0.5e0 + f*sx - g*cx;
				s = 0.5e0 - f*cx - g*sx;
				if (x < 0.0)
				{
					c = -c;
					s = -s;
				}
				havec = TRUE;
			}
			if (x < 0.0)
			{
				/*              We COULD do the following before the preceeding, and then
				 *              not put in a test in the preceeding for x .lt. 0, but
				 *              even though the results are mathematically identical, we
				 *              would have some cancellation above if we did so. */
				if (wantg)
					g = cx + sx - g;
				if (wantf)
					f = cx - sx - f;
			}
		}
	}
	xsave = x;
	switch (mode)
	{
		case 1: goto L_1;
		case 2: goto L_2;
		case 3: goto L_3;
		case 4: goto L_4;
	}
L_1:
	sfren1_v = c;
	return( sfren1_v );
L_2:
	sfren1_v = s;
	return( sfren1_v );
L_3:
	sfren1_v = f;
	return( sfren1_v );
L_4:
	sfren1_v = g;
	return( sfren1_v );
} /* end of function */