Geant4 Cross Reference |
1 /* 1 /* 2 # <<BEGIN-copyright>> 2 # <<BEGIN-copyright>> 3 # <<END-copyright>> 3 # <<END-copyright>> 4 */ 4 */ 5 5 6 #include <stdio.h> 6 #include <stdio.h> 7 #include <stdlib.h> 7 #include <stdlib.h> 8 #include <cmath> 8 #include <cmath> 9 #include <float.h> 9 #include <float.h> 10 10 11 #include "ptwXY.h" 11 #include "ptwXY.h" 12 12 13 #if defined __cplusplus 13 #if defined __cplusplus 14 #include <cmath> 14 #include <cmath> 15 #include "G4Log.hh" 15 #include "G4Log.hh" 16 namespace GIDI { 16 namespace GIDI { 17 using namespace GIDI; 17 using namespace GIDI; 18 #endif 18 #endif 19 19 20 static nfu_status ptwXY_createFromFunctionBise 20 static nfu_status ptwXY_createFromFunctionBisect( ptwXYPoints *ptwXY, double x1, double y1, double x2, double y2, ptwXY_createFromFunction_callback func, 21 void *argList, int level, int checkFor 21 void *argList, int level, int checkForRoots, double eps ); 22 static nfu_status ptwXY_createFromFunctionZero 22 static nfu_status ptwXY_createFromFunctionZeroCrossing( ptwXYPoints *ptwXY, double x1, double y1, double x2, double y2, 23 ptwXY_createFromFunction_callback func 23 ptwXY_createFromFunction_callback func, void *argList, double eps ); 24 static nfu_status ptwXY_applyFunction2( ptwXYP 24 static nfu_status ptwXY_applyFunction2( ptwXYPoints *ptwXY1, double y1, double y2, ptwXYPoint *p1, ptwXYPoint *p2, ptwXY_applyFunction_callback func, 25 void *argList, int level, int checkForRoot 25 void *argList, int level, int checkForRoots ); 26 static nfu_status ptwXY_applyFunctionZeroCross 26 static nfu_status ptwXY_applyFunctionZeroCrossing( ptwXYPoints *ptwXY1, double y1, double y2, ptwXYPoint *p1, ptwXYPoint *p2, 27 ptwXY_applyFunction_callback func, void *a 27 ptwXY_applyFunction_callback func, void *argList ); 28 /* 28 /* 29 ********************************************** 29 ************************************************************ 30 */ 30 */ 31 void ptwXY_update_biSectionMax( ptwXYPoints *p 31 void ptwXY_update_biSectionMax( ptwXYPoints *ptwXY1, double oldLength ) { 32 32 33 ptwXY1->biSectionMax = ptwXY1->biSectionMa 33 ptwXY1->biSectionMax = ptwXY1->biSectionMax - 1.442695 * G4Log( ptwXY1->length / oldLength ); /* 1.442695 = 1 / std::log( 2. ) */ 34 if( ptwXY1->biSectionMax < 0 ) ptwXY1->biS 34 if( ptwXY1->biSectionMax < 0 ) ptwXY1->biSectionMax = 0; 35 if( ptwXY1->biSectionMax > ptwXY_maxBiSect 35 if( ptwXY1->biSectionMax > ptwXY_maxBiSectionMax ) ptwXY1->biSectionMax = ptwXY_maxBiSectionMax; 36 } 36 } 37 /* 37 /* 38 ********************************************** 38 ************************************************************ 39 */ 39 */ 40 ptwXYPoints *ptwXY_createFromFunction( int n, 40 ptwXYPoints *ptwXY_createFromFunction( int n, double *xs, ptwXY_createFromFunction_callback func, void *argList, double accuracy, int checkForRoots, 41 int biSectionMax, nfu_status *status ) { 41 int biSectionMax, nfu_status *status ) { 42 42 43 int64_t i; 43 int64_t i; 44 double x1, y1, x2 = 0., y2, eps = ClosestA 44 double x1, y1, x2 = 0., y2, eps = ClosestAllowXFactor * DBL_EPSILON; 45 ptwXYPoints *ptwXY; 45 ptwXYPoints *ptwXY; 46 ptwXYPoint *p1, *p2; 46 ptwXYPoint *p1, *p2; 47 47 48 *status = nfu_Okay; 48 *status = nfu_Okay; 49 if( n < 2 ) { *status = nfu_tooFewPoints; 49 if( n < 2 ) { *status = nfu_tooFewPoints; return( NULL ); } 50 for( i = 1; i < n; i++ ) { 50 for( i = 1; i < n; i++ ) { 51 if( xs[i-1] >= xs[i] ) *status = nfu_X 51 if( xs[i-1] >= xs[i] ) *status = nfu_XNotAscending; 52 } 52 } 53 if( *status == nfu_XNotAscending ) return( 53 if( *status == nfu_XNotAscending ) return( NULL ); 54 54 55 x1 = xs[0]; 55 x1 = xs[0]; 56 if( ( *status = func( x1, &y1, argList ) ) 56 if( ( *status = func( x1, &y1, argList ) ) != nfu_Okay ) return( NULL ); 57 if( ( ptwXY = ptwXY_new( ptwXY_interpolati 57 if( ( ptwXY = ptwXY_new( ptwXY_interpolationLinLin, NULL, biSectionMax, accuracy, 500, 50, status, 0 ) ) == NULL ) return( NULL ); 58 for( i = 1; i < n; i++ ) { 58 for( i = 1; i < n; i++ ) { 59 if( ( *status = ptwXY_setValueAtX_over 59 if( ( *status = ptwXY_setValueAtX_overrideIfClose( ptwXY, x1, y1, eps, 0 ) ) != nfu_Okay ) goto err; 60 x2 = xs[i]; 60 x2 = xs[i]; 61 if( ( *status = func( x2, &y2, argList 61 if( ( *status = func( x2, &y2, argList ) ) != nfu_Okay ) goto err; 62 if( ( *status = ptwXY_createFromFuncti 62 if( ( *status = ptwXY_createFromFunctionBisect( ptwXY, x1, y1, x2, y2, func, argList, 0, checkForRoots, eps ) ) != nfu_Okay ) goto err; 63 x1 = x2; 63 x1 = x2; 64 y1 = y2; 64 y1 = y2; 65 } 65 } 66 if( ( *status = ptwXY_setValueAtX_override 66 if( ( *status = ptwXY_setValueAtX_overrideIfClose( ptwXY, x2, y2, eps, 1 ) ) != nfu_Okay ) goto err; 67 67 68 if( checkForRoots ) { 68 if( checkForRoots ) { 69 if( ( *status = ptwXY_simpleCoalescePo 69 if( ( *status = ptwXY_simpleCoalescePoints( ptwXY ) ) != nfu_Okay ) goto err; 70 for( i = ptwXY->length - 1, p2 = NULL; 70 for( i = ptwXY->length - 1, p2 = NULL; i >= 0; i--, p2 = p1 ) { /* Work backward so lower points are still valid if a new point is added. */ 71 p1 = &(ptwXY->points[i]); 71 p1 = &(ptwXY->points[i]); 72 if( p2 != NULL ) { 72 if( p2 != NULL ) { 73 if( ( p1->y * p2->y ) < 0. ) { 73 if( ( p1->y * p2->y ) < 0. ) { 74 if( ( *status = ptwXY_crea 74 if( ( *status = ptwXY_createFromFunctionZeroCrossing( ptwXY, p1->x, p1->y, p2->x, p2->y, func, argList, eps ) ) != nfu_Okay ) goto err; 75 } 75 } 76 } 76 } 77 } 77 } 78 } 78 } 79 79 80 return( ptwXY ); 80 return( ptwXY ); 81 81 82 err: 82 err: 83 ptwXY_free( ptwXY ); 83 ptwXY_free( ptwXY ); 84 return( NULL ); 84 return( NULL ); 85 } 85 } 86 /* 86 /* 87 ********************************************** 87 ************************************************************ 88 */ 88 */ 89 ptwXYPoints *ptwXY_createFromFunction2( ptwXPo 89 ptwXYPoints *ptwXY_createFromFunction2( ptwXPoints *xs, ptwXY_createFromFunction_callback func, void *argList, double accuracy, int checkForRoots, 90 int biSectionMax, nfu_status *status ) { 90 int biSectionMax, nfu_status *status ) { 91 91 92 return( ptwXY_createFromFunction( (int) xs 92 return( ptwXY_createFromFunction( (int) xs->length, xs->points, func, argList, accuracy, checkForRoots, biSectionMax, status ) ); 93 } 93 } 94 /* 94 /* 95 ********************************************** 95 ************************************************************ 96 */ 96 */ 97 static nfu_status ptwXY_createFromFunctionBise 97 static nfu_status ptwXY_createFromFunctionBisect( ptwXYPoints *ptwXY, double x1, double y1, double x2, double y2, ptwXY_createFromFunction_callback func, 98 void *argList, int level, int checkFor 98 void *argList, int level, int checkForRoots, double eps ) { 99 99 100 nfu_status status; 100 nfu_status status; 101 double x, y, f; 101 double x, y, f; 102 102 103 if( ( x2 - x1 ) < ClosestAllowXFactor * DB 103 if( ( x2 - x1 ) < ClosestAllowXFactor * DBL_EPSILON * ( std::fabs( x1 ) + std::fabs( x2 ) ) ) return( nfu_Okay ); 104 if( level >= ptwXY->biSectionMax ) return( 104 if( level >= ptwXY->biSectionMax ) return( nfu_Okay ); 105 x = 0.5 * ( x1 + x2 ); 105 x = 0.5 * ( x1 + x2 ); 106 if( ( status = ptwXY_interpolatePoint( ptw 106 if( ( status = ptwXY_interpolatePoint( ptwXY->interpolation, x, &y, x1, y1, x2, y2 ) ) != nfu_Okay ) return( status ); 107 if( ( status = func( x, &f, argList ) ) != 107 if( ( status = func( x, &f, argList ) ) != nfu_Okay ) return( status ); 108 if( std::fabs( f - y ) <= 0.8 * std::fabs( 108 if( std::fabs( f - y ) <= 0.8 * std::fabs( f * ptwXY->accuracy ) ) return( nfu_Okay ); 109 if( ( status = ptwXY_createFromFunctionBis 109 if( ( status = ptwXY_createFromFunctionBisect( ptwXY, x1, y1, x, f, func, argList, level + 1, checkForRoots, eps ) ) ) return( status ); 110 if( ( status = ptwXY_setValueAtX_overrideI 110 if( ( status = ptwXY_setValueAtX_overrideIfClose( ptwXY, x, f, eps, 0 ) ) != nfu_Okay ) return( status ); 111 return( ptwXY_createFromFunctionBisect( pt 111 return( ptwXY_createFromFunctionBisect( ptwXY, x, f, x2, y2, func, argList, level + 1, checkForRoots, eps ) ); 112 } 112 } 113 /* 113 /* 114 ********************************************** 114 ************************************************************ 115 */ 115 */ 116 static nfu_status ptwXY_createFromFunctionZero 116 static nfu_status ptwXY_createFromFunctionZeroCrossing( ptwXYPoints *ptwXY, double x1, double y1, double x2, double y2, 117 ptwXY_createFromFunction_callback func 117 ptwXY_createFromFunction_callback func, void *argList, double eps ) { 118 118 119 //For coverity #63077 119 //For coverity #63077 120 if ( y2 == y1 ) return ( nfu_badInput ); 120 if ( y2 == y1 ) return ( nfu_badInput ); 121 121 122 int i; 122 int i; 123 double x = 0, y = 0; << 123 double x, y; 124 nfu_status status; 124 nfu_status status; 125 125 126 for( i = 0; i < 16; i++ ) { 126 for( i = 0; i < 16; i++ ) { 127 if( y2 == y1 ) break; 127 if( y2 == y1 ) break; 128 x = ( y2 * x1 - y1 * x2 ) / ( y2 - y1 128 x = ( y2 * x1 - y1 * x2 ) / ( y2 - y1 ); 129 if( x <= x1 ) x = x1 + 0.1 * ( x2 - x1 129 if( x <= x1 ) x = x1 + 0.1 * ( x2 - x1 ); 130 if( x >= x2 ) x = x2 - 0.1 * ( x2 - x 130 if( x >= x2 ) x = x2 - 0.1 * ( x2 - x1 ); 131 if( ( status = func( x, &y, argList ) 131 if( ( status = func( x, &y, argList ) ) != nfu_Okay ) return( status ); 132 if( y == 0 ) break; 132 if( y == 0 ) break; 133 if( y1 * y < 0 ) { 133 if( y1 * y < 0 ) { 134 x2 = x; 134 x2 = x; 135 y2 = y; } 135 y2 = y; } 136 else { 136 else { 137 x1 = x; 137 x1 = x; 138 y1 = y; 138 y1 = y; 139 } 139 } 140 } 140 } 141 return( ptwXY_setValueAtX_overrideIfClose( 141 return( ptwXY_setValueAtX_overrideIfClose( ptwXY, x, 0., eps, 1 ) ); 142 } 142 } 143 /* 143 /* 144 ********************************************** 144 ************************************************************ 145 */ 145 */ 146 nfu_status ptwXY_applyFunction( ptwXYPoints *p 146 nfu_status ptwXY_applyFunction( ptwXYPoints *ptwXY1, ptwXY_applyFunction_callback func, void *argList, int checkForRoots ) { 147 147 148 int64_t i, originalLength = ptwXY1->length 148 int64_t i, originalLength = ptwXY1->length, notFirstPass = 0; 149 double y1, y2 = 0; 149 double y1, y2 = 0; 150 nfu_status status; 150 nfu_status status; 151 ptwXYPoint p1, p2; 151 ptwXYPoint p1, p2; 152 152 153 checkForRoots = checkForRoots && ptwXY1->b 153 checkForRoots = checkForRoots && ptwXY1->biSectionMax; 154 if( ptwXY1->status != nfu_Okay ) return( p 154 if( ptwXY1->status != nfu_Okay ) return( ptwXY1->status ); 155 if( ptwXY1->interpolation == ptwXY_interpo 155 if( ptwXY1->interpolation == ptwXY_interpolationOther ) return( nfu_otherInterpolation ); 156 if( ptwXY1->interpolation == ptwXY_interpo 156 if( ptwXY1->interpolation == ptwXY_interpolationFlat ) return( nfu_invalidInterpolation ); 157 if( ( status = ptwXY_simpleCoalescePoints( 157 if( ( status = ptwXY_simpleCoalescePoints( ptwXY1 ) ) != nfu_Okay ) return( status ); 158 for( i = originalLength - 1; i >= 0; i-- ) 158 for( i = originalLength - 1; i >= 0; i-- ) { 159 y1 = ptwXY1->points[i].y; 159 y1 = ptwXY1->points[i].y; 160 if( ( status = func( &(ptwXY1->points[ 160 if( ( status = func( &(ptwXY1->points[i]), argList ) ) != nfu_Okay ) return( status ); 161 p1 = ptwXY1->points[i]; 161 p1 = ptwXY1->points[i]; 162 if( notFirstPass ) { 162 if( notFirstPass ) { 163 if( ( status = ptwXY_applyFunction 163 if( ( status = ptwXY_applyFunction2( ptwXY1, y1, y2, &p1, &p2, func, argList, 0, checkForRoots ) ) != nfu_Okay ) return( status ); 164 } 164 } 165 notFirstPass = 1; 165 notFirstPass = 1; 166 p2 = p1; 166 p2 = p1; 167 y2 = y1; 167 y2 = y1; 168 } 168 } 169 ptwXY_update_biSectionMax( ptwXY1, (double 169 ptwXY_update_biSectionMax( ptwXY1, (double) originalLength ); 170 return( status ); 170 return( status ); 171 } 171 } 172 /* 172 /* 173 ********************************************** 173 ************************************************************ 174 */ 174 */ 175 static nfu_status ptwXY_applyFunction2( ptwXYP 175 static nfu_status ptwXY_applyFunction2( ptwXYPoints *ptwXY1, double y1, double y2, ptwXYPoint *p1, ptwXYPoint *p2, ptwXY_applyFunction_callback func, 176 void *argList, int level, int checkFor 176 void *argList, int level, int checkForRoots ) { 177 177 178 double y; 178 double y; 179 ptwXYPoint p; 179 ptwXYPoint p; 180 nfu_status status; 180 nfu_status status; 181 181 182 if( ( p2->x - p1->x ) < ClosestAllowXFacto 182 if( ( p2->x - p1->x ) < ClosestAllowXFactor * DBL_EPSILON * ( std::fabs( p1->x ) + std::fabs( p2->x ) ) ) return( nfu_Okay ); 183 if( level >= ptwXY1->biSectionMax ) goto c 183 if( level >= ptwXY1->biSectionMax ) goto checkForZeroCrossing; 184 p.x = 0.5 * ( p1->x + p2->x ); 184 p.x = 0.5 * ( p1->x + p2->x ); 185 if( ( status = ptwXY_interpolatePoint( ptw 185 if( ( status = ptwXY_interpolatePoint( ptwXY1->interpolation, p.x, &y, p1->x, y1, p2->x, y2 ) ) != nfu_Okay ) return( status ); 186 p.y = y; 186 p.y = y; 187 if( ( status = func( &p, argList ) ) != nf 187 if( ( status = func( &p, argList ) ) != nfu_Okay ) return( status ); 188 if( std::fabs( ( p.x - p1->x ) * ( p2->y - 188 if( std::fabs( ( p.x - p1->x ) * ( p2->y - p1->y ) + ( p2->x - p1->x ) * ( p1->y - p.y ) ) <= 0.8 * std::fabs( ( p2->x - p1->x ) * p.y * ptwXY1->accuracy ) ) 189 goto checkForZeroCrossing; 189 goto checkForZeroCrossing; 190 if( ( status = ptwXY_setValueAtX( ptwXY1, 190 if( ( status = ptwXY_setValueAtX( ptwXY1, p.x, p.y ) ) != nfu_Okay ) return( status ); 191 if( ( status = ptwXY_applyFunction2( ptwXY 191 if( ( status = ptwXY_applyFunction2( ptwXY1, y1, y, p1, &p, func, argList, level + 1, checkForRoots ) ) ) return( status ); 192 return( ptwXY_applyFunction2( ptwXY1, y, y 192 return( ptwXY_applyFunction2( ptwXY1, y, y2, &p, p2, func, argList, level + 1, checkForRoots ) ); 193 193 194 checkForZeroCrossing: 194 checkForZeroCrossing: 195 if( checkForRoots && ( ( p1->y * p2->y ) < 195 if( checkForRoots && ( ( p1->y * p2->y ) < 0. ) ) return( ptwXY_applyFunctionZeroCrossing( ptwXY1, y1, y2, p1, p2, func, argList ) ); 196 return( nfu_Okay ); 196 return( nfu_Okay ); 197 } 197 } 198 /* 198 /* 199 ********************************************** 199 ************************************************************ 200 */ 200 */ 201 static nfu_status ptwXY_applyFunctionZeroCross 201 static nfu_status ptwXY_applyFunctionZeroCrossing( ptwXYPoints *ptwXY1, double y1, double y2, ptwXYPoint *p1, ptwXYPoint *p2, 202 ptwXY_applyFunction_callback func, voi 202 ptwXY_applyFunction_callback func, void *argList ) { 203 203 204 int i; 204 int i; 205 double y, x1 = p1->x, x2 = p2->x, nY1 = p1 205 double y, x1 = p1->x, x2 = p2->x, nY1 = p1->y, nY2 = p2->y, refY = 0.5 * ( std::fabs( p1->y ) + std::fabs( p2->y ) ); 206 ptwXYPoint p; 206 ptwXYPoint p; 207 nfu_status status; 207 nfu_status status; 208 208 209 //For coverity #63074 209 //For coverity #63074 210 if ( nY2 == nY1 ) return ( nfu_badInput ); 210 if ( nY2 == nY1 ) return ( nfu_badInput ); 211 211 212 for( i = 0; i < 6; i++ ) { 212 for( i = 0; i < 6; i++ ) { 213 if( nY2 == nY1 ) break; 213 if( nY2 == nY1 ) break; 214 p.x = ( nY2 * x1 - nY1 * x2 ) / ( nY2 214 p.x = ( nY2 * x1 - nY1 * x2 ) / ( nY2 - nY1 ); 215 if( p.x <= x1 ) p.x = 0.5 * ( x1 + x2 215 if( p.x <= x1 ) p.x = 0.5 * ( x1 + x2 ); 216 if( p.x >= x2 ) p.x = 0.5 * ( x1 + x2 216 if( p.x >= x2 ) p.x = 0.5 * ( x1 + x2 ); 217 if( ( status = ptwXY_interpolatePoint( 217 if( ( status = ptwXY_interpolatePoint( ptwXY1->interpolation, p.x, &y, p1->x, y1, p2->x, y2 ) ) != nfu_Okay ) return( status ); 218 p.y = y; 218 p.y = y; 219 if( ( status = func( &p, argList ) ) ! 219 if( ( status = func( &p, argList ) ) != nfu_Okay ) return( status ); 220 if( p.y == 0 ) break; 220 if( p.y == 0 ) break; 221 if( 0.5 * refY < std::fabs( p.y ) ) br 221 if( 0.5 * refY < std::fabs( p.y ) ) break; 222 refY = std::fabs( p.y ); 222 refY = std::fabs( p.y ); 223 if( p1->y * p.y < 0 ) { 223 if( p1->y * p.y < 0 ) { 224 x2 = p.x; 224 x2 = p.x; 225 nY2 = p.y; } 225 nY2 = p.y; } 226 else { 226 else { 227 x1 = p.x; 227 x1 = p.x; 228 nY1 = p.y; 228 nY1 = p.y; 229 } 229 } 230 } 230 } 231 return( ptwXY_setValueAtX( ptwXY1, p.x, 0. 231 return( ptwXY_setValueAtX( ptwXY1, p.x, 0. ) ); 232 } 232 } 233 /* 233 /* 234 ********************************************** 234 ************************************************************ 235 */ 235 */ 236 ptwXYPoints *ptwXY_fromString( char const *str 236 ptwXYPoints *ptwXY_fromString( char const *str, ptwXY_interpolation interpolation, ptwXY_interpolationOtherInfo const *interpolationOtherInfo, 237 double biSectionMax, double accuracy, char 237 double biSectionMax, double accuracy, char **endCharacter, nfu_status *status ) { 238 238 239 int64_t numberConverted; 239 int64_t numberConverted; 240 double *doublePtr; 240 double *doublePtr; 241 ptwXYPoints *ptwXY = NULL; 241 ptwXYPoints *ptwXY = NULL; 242 242 243 if( ( *status = nfu_stringToListOfDoubles( 243 if( ( *status = nfu_stringToListOfDoubles( str, &numberConverted, &doublePtr, endCharacter ) ) != nfu_Okay ) return( NULL ); 244 *status = nfu_oddNumberOfValues; 244 *status = nfu_oddNumberOfValues; 245 if( ( numberConverted % 2 ) == 0 ) 245 if( ( numberConverted % 2 ) == 0 ) 246 ptwXY = ptwXY_create( interpolation, i 246 ptwXY = ptwXY_create( interpolation, interpolationOtherInfo, biSectionMax, accuracy, numberConverted, 10, numberConverted / 2, doublePtr, status, 0 ); 247 nfu_free( doublePtr ); 247 nfu_free( doublePtr ); 248 return( ptwXY ); 248 return( ptwXY ); 249 } 249 } 250 /* 250 /* 251 ********************************************** 251 ************************************************************ 252 */ 252 */ 253 void ptwXY_showInteralStructure( ptwXYPoints * 253 void ptwXY_showInteralStructure( ptwXYPoints *ptwXY, FILE *f, int printPointersAsNull ) { 254 254 255 int64_t i, n = ptwXY_getNonOverflowLength( 255 int64_t i, n = ptwXY_getNonOverflowLength( ptwXY ); 256 ptwXYPoint *point = ptwXY->points; 256 ptwXYPoint *point = ptwXY->points; 257 ptwXYOverflowPoint *overflowPoint; 257 ptwXYOverflowPoint *overflowPoint; 258 258 259 fprintf( f, "status = %d interpolation = 259 fprintf( f, "status = %d interpolation = %d length = %d allocatedSize = %d\n", 260 (int) ptwXY->status, (int) ptwXY->inte 260 (int) ptwXY->status, (int) ptwXY->interpolation, (int) ptwXY->length, (int) ptwXY->allocatedSize ); 261 fprintf( f, "userFlag = %d biSectionMax = 261 fprintf( f, "userFlag = %d biSectionMax = %.8e accuracy = %.2e minFractional_dx = %.6e\n", 262 ptwXY->userFlag, ptwXY->biSectionMax, 262 ptwXY->userFlag, ptwXY->biSectionMax, ptwXY->accuracy, ptwXY->minFractional_dx ); 263 fprintf( f, "interpolationString = %s\n", 263 fprintf( f, "interpolationString = %s\n", ptwXY->interpolationOtherInfo.interpolationString ); 264 fprintf( f, "getValueFunc is NULL = %d. ar 264 fprintf( f, "getValueFunc is NULL = %d. argList is NULL = %d.\n", 265 ptwXY->interpolationOtherInfo.getValue 265 ptwXY->interpolationOtherInfo.getValueFunc == NULL, ptwXY->interpolationOtherInfo.argList == NULL ); 266 fprintf( f, " overflowLength = %d overfl 266 fprintf( f, " overflowLength = %d overflowAllocatedSize = %d mallocFailedSize = %d\n", 267 (int) ptwXY->overflowLength, (int) ptw 267 (int) ptwXY->overflowLength, (int) ptwXY->overflowAllocatedSize, (int) ptwXY->mallocFailedSize ); 268 fprintf( f, " Points data, points = %20p\ 268 fprintf( f, " Points data, points = %20p\n", ( printPointersAsNull ? NULL : (void*)ptwXY->points ) ); 269 for( i = 0; i < n; i++, point++ ) fprintf 269 for( i = 0; i < n; i++, point++ ) fprintf( f, " %14.7e %14.7e\n", point->x, point->y ); 270 fprintf( f, " Overflow points data; %20p\ 270 fprintf( f, " Overflow points data; %20p\n", ( printPointersAsNull ? NULL : (void*)&(ptwXY->overflowHeader) ) ); 271 for( overflowPoint = ptwXY->overflowHeader 271 for( overflowPoint = ptwXY->overflowHeader.next; overflowPoint != &(ptwXY->overflowHeader); overflowPoint = overflowPoint->next ) { 272 fprintf( f, " %14.7e %14.7e %8d %20 272 fprintf( f, " %14.7e %14.7e %8d %20p %20p %20p\n", overflowPoint->point.x, overflowPoint->point.y, (int) overflowPoint->index, 273 (void*) ( printPointersAsNull ? NULL : ove 273 (void*) ( printPointersAsNull ? NULL : overflowPoint ), (void*) ( printPointersAsNull ? NULL : overflowPoint->prior ), 274 (void*) ( printPointersAsNull ? NULL : ove 274 (void*) ( printPointersAsNull ? NULL : overflowPoint->next ) ); 275 } 275 } 276 fprintf( f, " Points in order\n" ); 276 fprintf( f, " Points in order\n" ); 277 for( i = 0; i < ptwXY->length; i++ ) { 277 for( i = 0; i < ptwXY->length; i++ ) { 278 point = ptwXY_getPointAtIndex( ptwXY, 278 point = ptwXY_getPointAtIndex( ptwXY, i ); 279 fprintf( f, " %14.7e %14.7e\n", poi 279 fprintf( f, " %14.7e %14.7e\n", point->x, point->y ); 280 } 280 } 281 } 281 } 282 /* 282 /* 283 ********************************************** 283 ************************************************************ 284 */ 284 */ 285 void ptwXY_simpleWrite( ptwXYPoints *ptwXY, FI 285 void ptwXY_simpleWrite( ptwXYPoints *ptwXY, FILE *f, char *format ) { 286 286 287 int64_t i; 287 int64_t i; 288 ptwXYPoint *point; 288 ptwXYPoint *point; 289 289 290 for( i = 0; i < ptwXY->length; i++ ) { 290 for( i = 0; i < ptwXY->length; i++ ) { 291 point = ptwXY_getPointAtIndex( ptwXY, 291 point = ptwXY_getPointAtIndex( ptwXY, i ); 292 fprintf( f, format, point->x, point->y 292 fprintf( f, format, point->x, point->y ); 293 } 293 } 294 } 294 } 295 /* 295 /* 296 ********************************************** 296 ************************************************************ 297 */ 297 */ 298 void ptwXY_simplePrint( ptwXYPoints *ptwXY, ch 298 void ptwXY_simplePrint( ptwXYPoints *ptwXY, char *format ) { 299 299 300 ptwXY_simpleWrite( ptwXY, stdout, format ) 300 ptwXY_simpleWrite( ptwXY, stdout, format ); 301 } 301 } 302 302 303 #if defined __cplusplus 303 #if defined __cplusplus 304 } 304 } 305 #endif 305 #endif 306 306