/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:33:19 */ /*FOR_C Options SET: ftn=u io=c no=p op=aimnv pf=,p_ssfitc s=dbov str=l x=f - prototypes */ #include #include "fcrt.h" #include #include #include #include "p_ssfitc.h" /* program DRSSFITC *>> 2001-07-16 DRSSFITC Krogh Added exponent 0 to some constants. *>> 1996-07-11 DRSSFITC Krogh Special code for C conversion. *>> 1996-05-28 DRSSFITC Krogh Changed Fortran 90 code & changes for C. *>> 1994-10-19 DRSSFITC Krogh Changes to use M77CON *>> 1993-01-13 DRSSFITC C. L. Lawson, JPL *>> 1992-11-10 C. L. Lawson, JPL *>> 1992-11-04 C. L. Lawson, JPL *>> 1989-03-02 C. L. Lawson, JPL *>> 1988-04-01 C. L. Lawson, JPL * DRSSFITC.. Demo driver for SSFITC, Spline fit with constraints. * The problem has 24 data points and 10 constraints. * The spline is order 4 with 9 coefficients. * ------------------------------------------------------------------ *--S replaces "?": DR?SFITC, ?SFITC, ?SVAL, ?SVALA, ?PRPL, ?SDIF * ------------------------------------------------------------------ *++ Code for .C. is active */ long int k; #define MT (NCOEF+KORDER) #include /* PARAMETER translations */ #define KORDER 4 #define KPRINT 0 #define MXY (NDATA + 10) #define NCOEF 9 #define NDATA 24 #define NINFO 41 #define NWORK 843 /* end of PARAMETER translations */ int main( ) { char image[50]; long int i, info[NINFO], _i, _r; float bcoef[NCOEF], bdif[NCOEF*3], delx, rnorm, smax, smin, svalue[2-(0)+1], work[NWORK], x, yfit; static float sdi[MXY]; static float tknots[NCOEF + KORDER]={0.0e0,0.0e0,0.0e0,0.0e0,1.5e0, 2.5e0,3.3e0,4.0e0,4.7e0,6.0e0,6.0e0,6.0e0,6.0e0}; static char ccode[MXY + 1][5]={"10~a","10~a","10~a","10~a","10~a", "10~a","10~a","10~a","10~a","10~a","10~a","10~a","10~a","10~a", "10~a","10~a","10~a","10~a","10~a","10~a","10~a","10~a","10~a", "10~a","10=a","11>a","12>a","12>a","12>a","12a","10=a"," !"}; static float xi[MXY]={0.0e0,0.3e0,0.7e0,1.0e0,1.3e0,1.7e0,2.0e0, 2.3e0,2.5e0,2.6e0,2.8e0,2.9e0,3.0e0,3.1e0,3.2e0,3.5e0,3.7e0,4.0e0, 4.3e0,4.7e0,5.0e0,5.3e0,5.7e0,6.0e0,0.0e0,0.0e0,0.0e0,1.5e0,2.5e0, 3.5e0,4.5e0,6.0e0,6.0e0,6.0e0}; static float yi[MXY]={1.0e0,1.1e0,0.9e0,1.02e0,1.2e0,1.0e0,1.2e0, 1.4e0,1.76e0,2.0e0,2.4e0,2.6e0,3.0e0,3.4e0,3.7e0,4.3e0,4.45e0, 4.76e0,4.8e0,5.0e0,4.96e0,4.9e0,4.9e0,5.0e0,1.0e0,0.0e0,0.0e0, 0.0e0,0.0e0,0.0e0,0.0e0,0.0e0,0.0e0,5.0e0}; static long iset[3]={NINFO,NWORK,KPRINT}; static int _aini = 1; /* OFFSET Vectors w/subscript range: 1 to dimension */ float *const Bcoef = &bcoef[0] - 1; float *const Bdif = &bdif[0] - 1; long *const Info = &info[0] - 1; long *const Iset = &iset[0] - 1; float *const Sdi = &sdi[0] - 1; float *const Tknots = &tknots[0] - 1; float *const Work = &work[0] - 1; float *const Xi = &xi[0] - 1; float *const Yi = &yi[0] - 1; /* end of OFFSET VECTORS */ if( _aini ){ /* Do 1 TIME INITIALIZATIONS! */ Sdi[1] = -1.0e0; _aini = 0; } /*++ End */ /* ------------------------------------------------------------------ */ printf(" DRSSFITC.. Demo driver for SSFITC\n Least-squares polynomial spline fit to data with constraints.\n"); printf("\n I kind deriv relop active X Y\n"); for (i = 1; i <= MXY; i++) { printf(" %3ld %c %c %c %c%12.3f%10.3f\n", i, ccode[i - 1][0], ccode[i - 1][1], ccode[i - 1][2], ccode[i - 1][3], Xi[i], Yi[i]); } i = MXY + 1; printf(" %3ld %c %c %c %c\n", i, ccode[i - 1][0], ccode[i - 1][1], ccode[i - 1][2], ccode[i - 1][3]); printf("\n KORDER =%3d, NCOEF =%3d\n", KORDER, NCOEF); /*++ Code for ~.C. is inactive * print'('' TKNOTS() = '',4f10.5/(14x,4f10.5))', * * (TKNOTS(I), I = 1, MT) *++ Code for .C. is active */ printf("\n TKNOTS() = "); for (i = 1; i <= MT+3; i+=4){ for (k = i; k <= min( i+3, MT ); k++) printf("%10.5f", Tknots[k]); if (i + 3 < MT) printf("\n ");} printf("\n"); ssfitc( (byte(*)[5])ccode, xi, yi, sdi, KORDER, NCOEF, tknots, bcoef, &rnorm, iset, info, work ); /*++ End * * */ printf("\n After call to SSFITC:\n"); printf("\n IERR5 =%6ld, NEED1 =%7ld, NEED2 =%7ld\n M1 =%6ld, MFIT =%7ld, NS =%7ld\n RNORM =%12.5f\n", Info[1], Info[2], Info[3], Info[4], Info[5], Info[6], rnorm); /*++ Code for ~.C. is inactive * print'(/'' BCOEF() = '',4f10.5/(13x,4f10.5))', * * (BCOEF(I),I=1,NCOEF) *++ Code for .C. is active */ printf("\n BCOEF() = "); for (i = 1; i <= NCOEF+3; i+=4){ for (k = i; k <= min( i+3, NCOEF ); k++) printf("%10.5f", Bcoef[k]); if (i + 3 < NCOEF) printf("\n ");} printf("\n"); ssdif( KORDER, NCOEF, tknots, bcoef, 2, bdif ); /*++ End * */ smin = 0.0e0; smax = 0.0e0; delx = (Xi[NDATA] - Xi[1])/30.0e0; x = Xi[1]; for (i = 0; i <= 31; i++) { ssvala( KORDER, NCOEF, tknots, 2, bdif, x, svalue ); smin = fminf( fminf( smin, svalue[0] ), fminf( svalue[1], svalue[2] ) ); smax = fmaxf( fmaxf( smax, svalue[0] ), fmaxf( svalue[1], svalue[2] ) ); x += delx; } printf("\n X YFIT YFIT' YFIT'' \n"); x = Xi[1]; for (i = 0; i <= 31; i++) { ssvala( KORDER, NCOEF, tknots, 2, bdif, x, svalue ); strcpy( image, " " ); sprpl( svalue[0], '*', (byte*)image, 49, smin, smax, FALSE ); sprpl( svalue[1], '1', (byte*)image, 49, smin, smax, FALSE ); sprpl( svalue[2], '2', (byte*)image, 49, smin, smax, FALSE ); printf(" %6.3f%7.3f%7.3f%7.3f %49.49s\n", x, svalue[0], svalue[1], svalue[2], image); x += delx; } /* Compute and print residuals. */ printf("\n Residuals at the data points:\n\n I XI(I) YI(I) YFIT YFIT-YI(I) YFIT-YI(I)\n"); for (i = 1; i <= NDATA; i++) { yfit = ssval( KORDER, NCOEF, tknots, bcoef, Xi[i], 0 ); sprpl( yfit - Yi[i], '*', (byte*)image, 39, -0.18e0, 0.18e0, TRUE ); printf(" %4ld%8.3f%8.3f%8.3f%10.3f %39.39s\n", i, Xi[i], Yi[i], yfit, yfit - Yi[i], ntstr(image+0,39)); } exit(0); } /* end of function */