/*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 #include "fcrt.h" #include "stgrec.h" #include /* PARAMETER translations */ #define MXWRK 28 #define ZERO 0.0e0 /* end of PARAMETER translations */ void /*FUNCTION*/ stgrec( float x[], float y[], float z[], float dz[][2], long np, long triang[], long nt, long bdry[][4], long mb, float xylim[], long nx, long ny, float zfill, float *zvals, long mx, long my, long ncont, LOGICAL32 wantpd, float *dzvals) { #define ZVALS(I_,J_) (*(zvals+(I_)*(mx)+(J_))) #define DZVALS(I_,J_,K_) (*(dzvals+(I_)*(my)*(mx)+(J_)*(mx)+(K_))) LOGICAL32 first; long int i, ix, iy, mode; float delx, dely, dzout[2], svwrk0[MXWRK], svwrk1[MXWRK], xy[2]; /* OFFSET Vectors w/subscript range: 1 to dimension */ float *const Dzout = &dzout[0] - 1; float *const Svwrk0 = &svwrk0[0] - 1; float *const Svwrk1 = &svwrk1[0] - 1; float *const X = &x[0] - 1; float *const Xy = &xy[0] - 1; float *const Xylim = &xylim[0] - 1; float *const Y = &y[0] - 1; float *const Z = &z[0] - 1; /* end of OFFSET VECTORS */ /* Copyright (c) 1996 California Institute of Technology, Pasadena, CA. * ALL RIGHTS RESERVED. * Based on Government Sponsored Research NAS7-03001. *>> 1997-07-01 STGREC Krogh Reversed subscripts in B (CLL suggestion) *>> 1997-06-18 STGREC CLL Added arguments Bdry and MB. *>> 1996-02-02 STGREC CLL *>> 1996-01-11 STGREC CLL *>> 1995-11-03 STGREC CLL Added argument, NCONT. *>> 1995-09-26 STGREC CLL Editing for inclusion into MATH77. * This subr builds a rectangular array of Z values by interpo- * lation in a triangular grid. The subr will use either C0 or C1 * interpolation methods, depending on the setting of NCONT = 0 or 1. * Optionally this subr also builds a rectangular array of values of * first partial derivatives. * CLL, 1979 JUNE 12, MINOR CHANGES. * C.L.LAWSON, JPL, 1976 DEC 8. ADDED PARTIAL DERIVS 1977 MAR 3. * ------------------------------------------------------------------ * Subroutine arguments * * All args are of intent IN except ZVALS() and DZVALS() which are * of intent OUT. The contents of DZ() are only used when NCONT = 1. * * X(),Y(),Z() [floating,in] NP sets of (X,Y,Z) data. * * (DZ(I,J),I=1,2),J=1,NP) [floating,in] Partial derivs of Z w.r.t. * X and Y. The contents of this array are only used when NCONT = 1, * not when NCONT = 0. * * NP [integer,in] No. of data points. * * TRIANG(1:6*NT),NT [integer,in] TRIANG() contains pointers (integers) * describing the triangular grid. NT is the number of triangles. * Each triangle requires 6 pointers. * * Bdry(1:4, MB) [integer, in] Array containing pointers defining the * boundary of the (convex) triangular grid. * * MB [integer, in] Dimension parameter for Bdry() array. * * XYLIM(1:4) [floating,in] Contains XMIN, XMAX, YMIN, and YMAX * specifying boundaries of the rectangular grid. * * NX,NY [integer, in] NO. OF X AND Y GRID LINES IN RECTANGULAR GRID. * * ZFILL [floating, in] Selects action to be taken when a node of the * rectangular grid is outside the (convex) triangular grid. * If ZFILL = 0.0, extrapolated values will be provided by this * package. * If ZFILL .ne. 0.0, this subroutine will store ZFILL into ZVALS() * and DZVALS() for locations outside the triangular grid. * * ZVALS(,) [floating,out] MX BY MY ARRAY IN WHICH THE NX BY NY SET OF * INTERPOLATED VALUES WILL BE STORED. * * MX,MY [integer, in] DIMENSIONS OF THE ARRAYS ZVALS( , ) AND * DZVALS( , ,2) * * NCONT [integer,in] Order of continuity of interpolation method to be * used. * 0 selects a C0 method. This method uses linear interpolation over * each triangle. The function value is continuous across the whole * grid but generally the 1st partial derivatives are discontinuous * across boundaries between triangles. * 1 selects a C1 method. This method partitions each triangle into * 3 smaller triangles and constructs a bicubic patch over each of * the smaller triangles. The function value and 1st partial deriv- * atives are continuous across the whole grid but generally the 2nd * and higher order partial derivatives are discontinuous across * boundaries between the original triangles and between the small * triangles. * * WANTPD [logical, in] User sets = .true. if partial derivs are * desired, and .false. otherwise. * * DZVALS(,,) [floating,out] MX by NY by 2 array in which the * NX by NY by 2 set of partial derivatives will be stored if * WANTPD = .TRUE. * In the case of NCONT = 0, the interpolation function does not * generally have a partial derivative at points on the boundary * between two triangles. For such points this subr will return the * 1st partial derivatives associated with one of the adjoining * triangles. * ------------------------------------------------------------------ *--S replaces "?": ?TGREC, ?TGFI * ------------------------------------------------------------------ */ /* ------------------------------------------------------------------ */ delx = (Xylim[2] - Xylim[1])/(nx - 1); dely = (Xylim[4] - Xylim[3])/(ny - 1); for (i = 1; i <= MXWRK; i++) { Svwrk0[i] = ZERO; } /* XY(1) AND XY(2) WILL BE SET SUCCESSIVELY TO THE (X,Y) COORDS * OF EACH POINT OF THE RECTANGULAR GRID. * When NCONT .eq. 0, the values of the 1st derivatives are * constant within each triangle. Thus in the case of * NCONT .eq. 0 and WANTDZ being true, we only need to compute the * partial derivatives when NEWTRI is true. * */ Xy[2] = Xylim[3]; for (iy = 1; iy <= ny; iy++) { for (i = 1; i <= MXWRK; i++) { Svwrk1[i] = Svwrk0[i]; } first = TRUE; Xy[1] = Xylim[1]; for (ix = 1; ix <= nx; ix++) { stgfi( x, y, z, dz, triang, nt, bdry, mb, ncont, xy, &ZVALS(iy - 1,ix - 1), wantpd, dzout, &mode, svwrk1 ); if (first) { for (i = 1; i <= MXWRK; i++) { Svwrk0[i] = Svwrk1[i]; } first = FALSE; } if (mode == 0) { /* 0 means XY was found within the triangular grid. * */ if (wantpd) { DZVALS(0,iy - 1,ix - 1) = Dzout[1]; DZVALS(1,iy - 1,ix - 1) = Dzout[2]; } } else if (mode == 1) { /* 1 MEANS XY IS OUTSIDE THE TRIANGULAR GRID. * Will use extrapolated value if ZFILL .eq. zero. * Otherwise set results to the ZFILL * */ if (zfill == ZERO) { if (wantpd) { DZVALS(0,iy - 1,ix - 1) = Dzout[1]; DZVALS(1,iy - 1,ix - 1) = Dzout[2]; } } else { ZVALS(iy - 1,ix - 1) = zfill; if (wantpd) { DZVALS(0,iy - 1,ix - 1) = zfill; DZVALS(1,iy - 1,ix - 1) = zfill; } } } else { /* Mode = 2 MEANS CYCLING HAS HAPPENED IN the lookup * process. This should not happen. * ERROR MESSAGE WILL BE PRINTED FROM _TGFI. * */ ZVALS(iy - 1,ix - 1) = zfill; if (wantpd) { DZVALS(0,iy - 1,ix - 1) = zfill; DZVALS(1,iy - 1,ix - 1) = zfill; } } Xy[1] += delx; } Xy[2] += dely; } return; #undef ZVALS #undef DZVALS } /* end of function */