/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:30:10 */
/*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 "srlog.h"
#include <float.h>
#include <stdlib.h>
float /*FUNCTION*/ srlog(
float x)
{
	float r, srlog_v;
 
	/* Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
	 * ALL RIGHTS RESERVED.
	 * Based on Government Sponsored Research NAS7-03001.
	 *  File srlog contains user callable procedures srlog & srlog1, and
	 *               lower level procedure srlog2.
	 *  srlog(x)  computes x - 1 - ln(x).
	 *  srlog1(x) computes x - ln(1 + x).
	 *>> 1996-03-30 SRLOG  Krogh  Added external statements.
	 *>> 1994-11-09 CLL Edited to avoid ENTRY statement.
	 *>> 1994-10-20 SRLOG Krogh  Changes to use M77CON
	 *>> 1994-05-23 SRLOG WVS JPL Combine SRLOG and SRLOG1
	 *>> 1994-05-23 SRLOG WVS JPL Make SP and DP alike using CHGTYP
	 *>> 1993-05-06 SRLOG WVS JPL Conversion from NSWC to Math 77
	 * ----------------------------------------------------------------------
	 *--S replaces "?": ?RLOG, ?RLOG1, ?RLOG2
	 *     ==================================================================
	 *             EVALUATION OF THE FUNCTION X - 1 - LN(X)
	 *     ------------------------------------------------------------------ */
	/*     ------------------------------------------------------------------ */
	r = (x - 0.5e0) - 0.5e0;
	if (x < 0.61e0 || x > 1.57e0)
	{
		srlog_v = r - logf( x );
	}
	else
	{
		srlog_v = srlog2( r );
	}
	return( srlog_v );
} /* end of function */
/*     ================================================================== */
float /*FUNCTION*/ srlog1(
float x)
{
	float r, srlog1_v;
 
	/*             EVALUATION OF THE FUNCTION X - LN(1 + X)
	 *     ------------------------------------------------------------------ */
	/*     ------------------------------------------------------------------ */
	if (x < -0.39e0 || x > 0.57e0)
	{
		r = (x + 0.5e0) + 0.5e0;
		srlog1_v = x - logf( r );
	}
	else
	{
		srlog1_v = srlog2( x );
	}
	return( srlog1_v );
} /* end of function */
/*     ================================================================== */
		/* PARAMETER translations */
#define	C1	(1.0e0/3.0e0)
#define	C2	(1.0e0/5.0e0)
#define	C3	(1.0e0/7.0e0)
#define	C4	(1.0e0/9.0e0)
#define	C5	(1.0e0/11.0e0)
		/* end of PARAMETER translations */
 
float /*FUNCTION*/ srlog2(
float rin)
{
	long int _l0;
	float r, srlog2_v, t, u, up2, w, z;
	static float a = .566749439387323789126387112411845e-01;
	static float b = .456512608815524058941143273395059e-01;
	static float p0 = .7692307692307692307680e-01;
	static float p1 = -.1505958055914600184836e00;
	static float p2 = .9302355725278521726994e-01;
	static float p3 = -.1787900022182327735804e-01;
	static float q1 = -.2824412139355646910683e01;
	static float q2 = .2892424216041495392509e01;
	static float q3 = -.1263560605948009364422e01;
	static float q4 = .1966769435894561313526e00;
	static float r0 = .333333333333333e00;
	static float r1 = -.224696413112536e00;
	static float r2 = .620886815375787e-02;
	static float s1 = -.127408923933623e01;
	static float s2 = .354508718369557e00;
	static float round = -1.0e0;
 
	/*            Complete computation started by SRLOG or SRLOG1.
	 *     ------------------------------------------------------------------ */
	/* ------------------------
	 *     CI = 1/(2I + 1)
	 * ------------------------ */
	/* ------------------------
	 *     A = SRLOG (0.7)
	 *     B = SRLOG (4/3)
	 * ------------------------ */
	/* ------------------------ */
	/* ------------------------ */
	/* ------------------------ */
	/*     ------------------------------------------------------------------
	 *
	 *                 ARGUMENT REDUCTION
	 * */
	r = rin;
	if (r < -0.18e0)
	{
		u = (10.0e0*r + 3.0e0)/7.0e0;
		up2 = u + 2.0e0;
		w = a - u*0.3e0;
	}
	else if (r > 0.18e0)
	{
		t = 0.75e0*r;
		u = t - 0.25e0;
		up2 = t + 1.75e0;
		w = b + u/3.0e0;
	}
	else
	{
		u = r;
		up2 = u + 2.0e0;
		w = 0.0e0;
	}
	if (round < 0.0e0)
		round = FLT_EPSILON;
 
	/*                  SERIES EXPANSION
	 * */
	r = u/up2;
	t = r*r;
	if (round > 5.0e-14)
	{
		z = ((r2*t + r1)*t + r0)/((s2*t + s1)*t + 1.0e0);
	}
	else
	{
 
		/*        Z IS (AT FIRST) A MINIMAX APPROXIMATION OF THE SERIES
		 *
		 *               C6 + C7*R**2 + C8*R**4 + ...
		 *
		 *        FOR THE INTERVAL (0.0, 0.375). THE APPROXIMATION IS ACCURATE
		 *        TO WITHIN 1.6 UNITS OF THE 21-ST SIGNIFICANT DIGIT.
		 * */
		z = (((p3*t + p2)*t + p1)*t + p0)/((((q4*t + q3)*t + q2)*t +
		 q1)*t + 1.0e0);
 
		z = ((((z*t + C5)*t + C4)*t + C3)*t + C2)*t + C1;
	}
	srlog2_v = r*(u - 2.0e0*t*z) + w;
	return( srlog2_v );
} /* end of function */