/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:31:48 */
/*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 "ssfind.h"
#include <stdlib.h>
void /*FUNCTION*/ ssfind(
float xt[],
long ix1,
long ix2,
float x,
long *lefti,
long *mode)
{
long int ihi, ilo, middle, step;
/* Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
* ALL RIGHTS RESERVED.
* Based on Government Sponsored Research NAS7-03001.
*>> 1996-03-30 SSFIND Krogh Removed Fortran 90 comment.
*>> 1994-11-11 CLL Replaced Fortran 90 DO WHILE's with GOTO's for F77.
*>> 1994-10-20 SSFIND Krogh Changes to use M77CON
*>> 1992-11-12 SSFIND C. L. Lawson, JPL
*>> 1992-10-23 C. L. Lawson, JPL
*>> 1988-03-16 C. L. Lawson, JPL
*
* Require IX1 < IX2.
* Require the values in XT() indexed from IX1 to IX2 to be
* nondecreasing.
* Require XT(IX1) < XT(IX1+1).
* Require XT(IX2-1) < XT(IX2).
* Let A = XT(IX1) and B = XT(IX2).
* The closed interval, [A, B], is regarded as being partitioned
* partitioned into IX2-IX1 disjoint subintervals, with all but the
* last being half-open, and the last one being closed:
* [XT(I), XT(I+1)), I = IX1, ..., IX2-2
* [XT(IX2-1), XT(IX2)]
* Some of these intervals may have zero length, but not the first
* and last ones.
* This subroutine identifies the location of X with respect to these
* subintervals by setting LEFTI and MODE as follows:
* If X is contained in one of these segments of nonzero length,
* set MODE = 0 and set LEFTI to be the index of the left end of the
* segment containing X. Thus LEFTI will satisfy
* IX1 .le. LEFTI .le. IX2-1.
* If X < A, set MODE = -1 and LEFTI = IX1.
* If X > B, set MODE = +1 and LEFTI = IX2-1.
* Issue an error message and stop if
* X < XT(IX1+1) and XT(IX1) .ge. XT(IX1+1) or if
* X .ge. XT(IX2-1) and XT(IX2-1) .ge. XT(IX2)
* ------------------------------------------------------------------
* Method
*
* Saves the value of LEFTI returned on previous call in ILO.
* Starts by checking to see if X is in this same segment.
* If so, finishes quickly. We assume this will be a frequently
* occurring case.
* If not, searches either to the left or right from this previous
* segment, as appropriate. During this search the stride is doubled
* until a bracketing value is found. Then we finish with a binary
* search between the bracketing points.
* ------------------------------------------------------------------
* Based on subroutine INTERV on pp. 92-93 of A PRACTICAL GUIDE TO
* SPLINES by Carl De Boor, Springer-Verlag, 1978.
* Current version by C. L. Lawson, JPL, March 1988.
* ------------------------------------------------------------------
* XT() [in] This subr will access XT(i) for i = IX1, ..., IX2. These
* values must be nondecreasing. Repeated values are permitted
* except for the first two and last two.
* IX1, IX2 [in] Indices specifying the portion of the array XT() to
* be considered by this subr. Require IX1 < IX2.
* X [in] Value to be looked up.
* LEFTI [inout] On entry LEFTI must contain an integer value. If it
* in the range [IX1, IX2-1] the search will start with this index.
* Otherwise the search will start with IX1 or IX2-1.
* On return LEFTI will be in [IX1, IX2-1], and identifies
* the interval (of nonzero length) from X(LEFTI) to
* X(LEFTI+1) as the reference subinterval for X.
* LEFTI = IX1 means X < XT(IX1+1)
* IX1 < LEFTI < IX2-1 means XT(LEFTI) .le. X .lt. XT(LEFTI+1)
* LEFTI = IX2-1 means XT(IX2-1) .le. X
* MODE [out] Set to
* -1 if X < XT(IX1)
* 0 if XT(IX1) .le. X .le. XT(IX2)
* +1 if XT(IX2) < X
* ------------------------------------------------------------------
*--S replaces "?": ?SFIND, ?ERV1
* Both use IERM1
* ------------------------------------------------------------------ */
/* ------------------------------------------------------------------ */
ilo = max( ix1, min( *lefti, ix2 - 1 ) );
if (x >= xt[ilo-(1)])
{
if (x < xt[ilo + 1-(1)])
{
*mode = 0;
*lefti = ilo;
return;
}
else
{
/* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
* procedure( SEARCH TO RIGHT ) */
step = 1;
ihi = ilo + 1;
/* do while(.true.) ! Using Fortran 90 "DO" syntax. */
L_5:
;
if (ihi == ix2)
{
/* Here X .ge. XT(IX2) */
if (x == xt[ix2-(1)])
{
*mode = 0;
}
else
{
*mode = 1;
}
ilo = ix2 - 1;
*lefti = ilo;
if (xt[ilo-(1)] >= xt[ix2-(1)])
{
ierm1( "SSFIND", 2, 2, "Require T(IX2-1) < T(IX2)"
, "IX2", ix2, ',' );
serv1( "T(IX2-1)", xt[ix2 - 1-(1)], ',' );
serv1( "T(IX2)", xt[ix2-(1)], '.' );
}
return;
}
ilo = ihi;
ihi = min( ihi + step, ix2 );
if (x < xt[ihi-(1)])
goto L_10;
step *= 2.0e0;
goto L_5;
/* end do !while */
L_10:
;
/* end proc !( SEARCH TO RIGHT )
* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
}
}
else
{
/* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
* procedure( SEARCH TO LEFT ) */
step = 1;
/* do while(.true.) */
L_15:
;
if (ilo == ix1)
{
*mode = -1;
*lefti = ix1;
if (xt[ix1-(1)] >= xt[ix1 + 1-(1)])
{
ierm1( "SSFIND", 1, 2, "Require T(IX1) < T(IX1+1)"
, "IX1", ix1, ',' );
serv1( "T(IX1)", xt[ix1-(1)], ',' );
serv1( "T(IX1+1)", xt[ix1 + 1-(1)], '.' );
}
return;
}
ihi = ilo;
ilo = max( ilo - step, ix1 );
if (x >= xt[ilo-(1)])
goto L_20;
step *= 2.0e0;
goto L_15;
/* end do !while */
L_20:
;
/* end proc !( SEARCH TO LEFT )
* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
}
/* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
* procedure( BINARY SEARCH ) */
/* Here XT(ILO) .le. X .lt. XT(IHI)
*
* do while(.true.) */
L_25:
;
middle = (ilo + ihi)/2;
if (middle == ilo)
goto L_30;
if (x >= xt[middle-(1)])
{
ilo = middle;
}
else
{
ihi = middle;
}
goto L_25;
/* end do !while */
L_30:
;
*mode = 0;
*lefti = ilo;
/* end proc !( BINARY SEARCH )
* -- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
return;
} /* end of function */