PLplot 5.9.6
|
00001 /* $Id$ 00002 * 00003 * Functions to shade regions on the basis of value. 00004 * Can be used to shade contour plots or alone. 00005 * Copyright 1993 Wesley Ebisuzaki 00006 * 00007 * Copyright (C) 2004 Alan W. Irwin 00008 * 00009 * This file is part of PLplot. 00010 * 00011 * PLplot is free software; you can redistribute it and/or modify 00012 * it under the terms of the GNU General Library Public License as published 00013 * by the Free Software Foundation; either version 2 of the License, or 00014 * (at your option) any later version. 00015 * 00016 * PLplot is distributed in the hope that it will be useful, 00017 * but WITHOUT ANY WARRANTY; without even the implied warranty of 00018 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00019 * GNU Library General Public License for more details. 00020 * 00021 * You should have received a copy of the GNU Library General Public License 00022 * along with PLplot; if not, write to the Free Software 00023 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 00024 * 00025 */ 00026 00027 /*----------------------------------------------------------------------*\ 00028 * Call syntax for plshade(): 00029 * 00030 * void plshade(PLFLT *a, PLINT nx, PLINT ny, char *defined, 00031 * PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00032 * PLFLT shade_min, PLFLT shade_max, 00033 * PLINT sh_color, PLINT sh_width, PLINT min_color, PLINT min_width, 00034 * PLINT max_color, PLINT max_width, void (*fill)(), PLINT 00035 * rectangular, ...) 00036 * 00037 * arguments: 00038 * 00039 * PLFLT &(a[0][0]) 00040 * 00041 * Contains array to be plotted. The array must have been declared as 00042 * PLFLT a[nx][ny]. See following note on fortran-style arrays. 00043 * 00044 * PLINT nx, ny 00045 * 00046 * Dimension of array "a". 00047 * 00048 * char &(defined[0][0]) 00049 * 00050 * Contains array of flags, 1 = data is valid, 0 = data is not valid. 00051 * This array determines which sections of the data is to be plotted. 00052 * This argument can be NULL if all the values are valid. Must have been 00053 * declared as char defined[nx][ny]. 00054 * 00055 * PLFLT xmin, xmax, ymin, ymax 00056 * 00057 * Defines the "grid" coordinates. The data a[0][0] has a position of 00058 * (xmin,ymin). 00059 * 00060 * void (*mapform)() 00061 * 00062 * Transformation from `grid' coordinates to world coordinates. This 00063 * pointer to a function can be NULL in which case the grid coordinates 00064 * are the same as the world coordinates. 00065 * 00066 * PLFLT shade_min, shade_max 00067 * 00068 * Defines the interval to be shaded. If shade_max <= shade_min, plshade 00069 * does nothing. 00070 * 00071 * PLINT sh_cmap, PLFLT sh_color, PLINT sh_width 00072 * 00073 * Defines color map, color map index, and width used by the fill pattern. 00074 * 00075 * PLINT min_color, min_width, max_color, max_width 00076 * 00077 * Defines pen color, width used by the boundary of shaded region. The min 00078 * values are used for the shade_min boundary, and the max values are used 00079 * on the shade_max boundary. Set color and width to zero for no plotted 00080 * boundaries. 00081 * 00082 * void (*fill)() 00083 * 00084 * Routine used to fill the region. Use plfill. Future version of plplot 00085 * may have other fill routines. 00086 * 00087 * PLINT rectangular 00088 * 00089 * Flag. Set to 1 if rectangles map to rectangles after (*mapform)() else 00090 * set to zero. If rectangular is set to 1, plshade tries to save time by 00091 * filling large rectangles. This optimization fails if (*mapform)() 00092 * distorts the shape of rectangles. For example a plot in polor 00093 * coordinates has to have rectangular set to zero. 00094 * 00095 * Example mapform's: 00096 * 00097 * Grid to world coordinate transformation. 00098 * This example goes from a r-theta to x-y for a polar plot. 00099 * 00100 * void mapform(PLINT n, PLFLT *x, PLFLT *y) { 00101 * int i; 00102 * double r, theta; 00103 * for (i = 0; i < n; i++) { 00104 * r = x[i]; 00105 * theta = y[i]; 00106 * x[i] = r*cos(theta); 00107 * y[i] = r*sin(theta); 00108 * } 00109 * } 00110 * 00111 * Grid was in cm, convert to world coordinates in inches. 00112 * Expands in x direction. 00113 * 00114 * void mapform(PLINT n, PLFLT *x, PLFLT *y) { 00115 * int i; 00116 * for (i = 0; i < n; i++) { 00117 * x[i] = (1.0 / 2.5) * x[i]; 00118 * y[i] = (1.0 / 2.5) * y[i]; 00119 * } 00120 * } 00121 * 00122 \*----------------------------------------------------------------------*/ 00123 00124 #include "plplotP.h" 00125 #include <float.h> 00126 00127 #define MISSING_MIN_DEF (PLFLT) 1.0 00128 #define MISSING_MAX_DEF (PLFLT) -1.0 00129 00130 00131 #define NEG 1 00132 #define POS 8 00133 #define OK 0 00134 #define UNDEF 64 00135 #define NUMBER_BISECTIONS 10 00136 00137 #define linear( val1, val2, level ) ( ( level - val1 ) / ( val2 - val1 ) ) 00138 00139 /* Global variables */ 00140 00141 static PLFLT sh_max, sh_min; 00142 static int min_points, max_points, n_point; 00143 static int min_pts[4], max_pts[4]; 00144 static PLINT pen_col_min, pen_col_max; 00145 static PLINT pen_wd_min, pen_wd_max; 00146 static PLFLT int_val; 00147 00148 /* Function prototypes */ 00149 00150 static void 00151 set_cond( register int *cond, register PLFLT *a, register PLINT n ); 00152 00153 static int 00154 find_interval( PLFLT a0, PLFLT a1, PLINT c0, PLINT c1, PLFLT *x ); 00155 00156 static void 00157 selected_polygon( void ( *fill )( PLINT, PLFLT *, PLFLT * ), 00158 PLINT ( *defined )( PLFLT, PLFLT ), 00159 PLFLT *x, PLFLT *y, PLINT v1, PLINT v2, PLINT v3, PLINT v4 ); 00160 00161 static void 00162 exfill( void ( *fill )( PLINT, PLFLT *, PLFLT * ), 00163 PLINT ( *defined )( PLFLT, PLFLT ), 00164 int n, PLFLT *x, PLFLT *y ); 00165 00166 static void 00167 big_recl( int *cond_code, register int ny, int dx, int dy, 00168 int *ix, int *iy ); 00169 00170 static void 00171 draw_boundary( PLINT slope, PLFLT *x, PLFLT *y ); 00172 00173 static PLINT 00174 plctest( PLFLT *x, PLFLT level ); 00175 00176 static PLINT 00177 plctestez( PLFLT *a, PLINT nx, PLINT ny, PLINT ix, 00178 PLINT iy, PLFLT level ); 00179 00180 static void 00181 plshade_int( PLFLT ( *f2eval )( PLINT, PLINT, PLPointer ), 00182 PLPointer f2eval_data, 00183 PLFLT ( *c2eval )( PLINT, PLINT, PLPointer ), 00184 PLPointer c2eval_data, 00185 PLINT ( *defined )( PLFLT, PLFLT ), 00186 PLFLT missing_min, PLFLT missing_max, 00187 PLINT nx, PLINT ny, 00188 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00189 PLFLT shade_min, PLFLT shade_max, 00190 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00191 PLINT min_color, PLINT min_width, 00192 PLINT max_color, PLINT max_width, 00193 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00194 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00195 PLPointer pltr_data ); 00196 00197 /*----------------------------------------------------------------------*\ 00198 * plshades() 00199 * 00200 * Shade regions via a series of calls to plshade. 00201 * All arguments are the same as plshade except the following: 00202 * clevel is a pointer to an array of values representing 00203 * the shade edge values, nlevel-1 is 00204 * the number of different shades, (nlevel is the number of shade edges), 00205 * fill_width is the pattern fill width, and cont_color and cont_width 00206 * are the color and width of the contour drawn at each shade edge. 00207 * (if cont_color <= 0 or cont_width <=0, no such contours are drawn). 00208 \*----------------------------------------------------------------------*/ 00209 00210 void c_plshades( PLFLT **a, PLINT nx, PLINT ny, PLINT ( *defined )( PLFLT, PLFLT ), 00211 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00212 PLFLT *clevel, PLINT nlevel, PLINT fill_width, 00213 PLINT cont_color, PLINT cont_width, 00214 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00215 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00216 PLPointer pltr_data ) 00217 { 00218 plfshades( plf2ops_c(), a, nx, ny, defined, 00219 xmin, xmax, ymin, ymax, 00220 clevel, nlevel, fill_width, 00221 cont_color, cont_width, 00222 fill, rectangular, 00223 pltr, pltr_data ); 00224 } 00225 00226 /*----------------------------------------------------------------------*\ 00227 * plfshades() 00228 * 00229 * Shade regions via a series of calls to plfshade1. 00230 * All arguments are the same as plfshade1 except the following: 00231 * clevel is a pointer to an array of values representing 00232 * the shade edge values, nlevel-1 is 00233 * the number of different shades, (nlevel is the number of shade edges), 00234 * fill_width is the pattern fill width, and cont_color and cont_width 00235 * are the color and width of the contour drawn at each shade edge. 00236 * (if cont_color <= 0 or cont_width <=0, no such contours are drawn). 00237 \*----------------------------------------------------------------------*/ 00238 00239 void 00240 plfshades( PLF2OPS zops, PLPointer zp, PLINT nx, PLINT ny, 00241 PLINT ( *defined )( PLFLT, PLFLT ), 00242 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00243 PLFLT *clevel, PLINT nlevel, PLINT fill_width, 00244 PLINT cont_color, PLINT cont_width, 00245 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00246 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00247 PLPointer pltr_data ) 00248 { 00249 PLFLT shade_min, shade_max, shade_color; 00250 PLINT i, init_color, init_width; 00251 00252 for ( i = 0; i < nlevel - 1; i++ ) 00253 { 00254 shade_min = clevel[i]; 00255 shade_max = clevel[i + 1]; 00256 shade_color = i / (PLFLT) ( nlevel - 2 ); 00257 /* The constants in order mean 00258 * (1) color map1, 00259 * (0, 0, 0, 0) all edge effects will be done with plcont rather 00260 * than the normal plshade drawing which gets partially blocked 00261 * when sequential shading is done as in the present case */ 00262 00263 plfshade1( zops, zp, nx, ny, defined, xmin, xmax, ymin, ymax, 00264 shade_min, shade_max, 00265 1, shade_color, fill_width, 00266 0, 0, 0, 0, 00267 fill, rectangular, pltr, pltr_data ); 00268 } 00269 if ( cont_color > 0 && cont_width > 0 ) 00270 { 00271 init_color = plsc->icol0; 00272 init_width = plsc->width; 00273 plcol0( cont_color ); 00274 plwid( cont_width ); 00275 if ( pltr ) 00276 { 00277 plfcont( zops->f2eval, zp, nx, ny, 1, nx, 1, ny, clevel, nlevel, pltr, pltr_data ); 00278 } 00279 else 00280 { 00281 /* For this case use the same interpretation that occurs internally 00282 * for plshade. That is set up x and y grids that map from the 00283 * index ranges to xmin, xmax, ymin, ymax, and use those grids 00284 * for the plcont call. 00285 */ 00286 PLcGrid cgrid1; 00287 PLFLT *x, *y; 00288 cgrid1.nx = nx; 00289 cgrid1.ny = ny; 00290 x = (PLFLT *) malloc( nx * sizeof ( PLFLT ) ); 00291 if ( x == NULL ) 00292 plexit( "plfshades: Out of memory for x" ); 00293 cgrid1.xg = x; 00294 for ( i = 0; i < nx; i++ ) 00295 cgrid1.xg[i] = xmin + ( xmax - xmin ) * (float) i / (float) ( nx - 1 ); 00296 y = (PLFLT *) malloc( ny * sizeof ( PLFLT ) ); 00297 if ( y == NULL ) 00298 plexit( "plfshades: Out of memory for y" ); 00299 cgrid1.yg = y; 00300 for ( i = 0; i < ny; i++ ) 00301 cgrid1.yg[i] = ymin + ( ymax - ymin ) * (float) i / (float) ( ny - 1 ); 00302 plfcont( zops->f2eval, zp, nx, ny, 1, nx, 1, ny, clevel, nlevel, 00303 pltr1, (void *) &cgrid1 ); 00304 free( x ); 00305 free( y ); 00306 } 00307 plcol0( init_color ); 00308 plwid( init_width ); 00309 } 00310 } 00311 00312 /*----------------------------------------------------------------------*\ 00313 * plshade() 00314 * 00315 * Shade region. 00316 * This interface to plfshade() assumes the 2d function array is passed 00317 * via a (PLFLT **), and is column-dominant (normal C ordering). 00318 \*----------------------------------------------------------------------*/ 00319 00320 void c_plshade( PLFLT **a, PLINT nx, PLINT ny, PLINT ( *defined )( PLFLT, PLFLT ), 00321 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00322 PLFLT shade_min, PLFLT shade_max, 00323 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00324 PLINT min_color, PLINT min_width, 00325 PLINT max_color, PLINT max_width, 00326 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00327 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00328 PLPointer pltr_data ) 00329 { 00330 plshade_int( plf2eval1, ( PLPointer ) a, 00331 NULL, NULL, 00332 /* plc2eval, (PLPointer) &cgrid,*/ 00333 defined, MISSING_MIN_DEF, MISSING_MAX_DEF, nx, ny, xmin, 00334 xmax, ymin, ymax, shade_min, shade_max, 00335 sh_cmap, sh_color, sh_width, 00336 min_color, min_width, max_color, max_width, 00337 fill, rectangular, pltr, pltr_data ); 00338 } 00339 00340 /*----------------------------------------------------------------------*\ 00341 * plshade1() 00342 * 00343 * Shade region. 00344 * This interface to plfshade() assumes the 2d function array is passed 00345 * via a (PLFLT *), and is column-dominant (normal C ordering). 00346 \*----------------------------------------------------------------------*/ 00347 00348 void c_plshade1( PLFLT *a, PLINT nx, PLINT ny, PLINT ( *defined )( PLFLT, PLFLT ), 00349 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00350 PLFLT shade_min, PLFLT shade_max, 00351 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00352 PLINT min_color, PLINT min_width, 00353 PLINT max_color, PLINT max_width, 00354 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00355 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00356 PLPointer pltr_data ) 00357 { 00358 PLfGrid grid; 00359 00360 grid.f = a; 00361 grid.nx = nx; 00362 grid.ny = ny; 00363 00364 plshade_int( plf2eval, ( PLPointer ) & grid, 00365 NULL, NULL, 00366 /* plc2eval, (PLPointer) &cgrid,*/ 00367 defined, MISSING_MIN_DEF, MISSING_MAX_DEF, nx, ny, xmin, 00368 xmax, ymin, ymax, shade_min, shade_max, 00369 sh_cmap, sh_color, sh_width, 00370 min_color, min_width, max_color, max_width, 00371 fill, rectangular, pltr, pltr_data ); 00372 } 00373 00374 /*----------------------------------------------------------------------*\ 00375 * plfshade() 00376 * 00377 * Shade region. 00378 * Array values are determined by the input function and the passed data. 00379 \*----------------------------------------------------------------------*/ 00380 00381 void 00382 plfshade( PLFLT ( *f2eval )( PLINT, PLINT, PLPointer ), 00383 PLPointer f2eval_data, 00384 PLFLT ( *c2eval )( PLINT, PLINT, PLPointer ), 00385 PLPointer c2eval_data, 00386 PLINT nx, PLINT ny, 00387 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00388 PLFLT shade_min, PLFLT shade_max, 00389 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00390 PLINT min_color, PLINT min_width, 00391 PLINT max_color, PLINT max_width, 00392 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00393 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00394 PLPointer pltr_data ) 00395 { 00396 plshade_int( f2eval, f2eval_data, c2eval, c2eval_data, 00397 NULL, MISSING_MIN_DEF, MISSING_MAX_DEF, 00398 nx, ny, xmin, xmax, ymin, ymax, 00399 shade_min, shade_max, sh_cmap, sh_color, sh_width, 00400 min_color, min_width, max_color, max_width, 00401 fill, rectangular, pltr, pltr_data ); 00402 } 00403 00404 /*----------------------------------------------------------------------*\ 00405 * plfshade1() 00406 * 00407 * Shade region. 00408 * 00409 * This function is a plf2ops variant of c_plfshade and c_plfshade1. It 00410 * differs from plfshade in that it supports a "defined" callback (like 00411 * c_plshade and c_plfshade1) rather than a "defined mask" (like plfshade 00412 * even though it is not yet implemented). 00413 \*----------------------------------------------------------------------*/ 00414 00415 void 00416 plfshade1( PLF2OPS zops, PLPointer zp, PLINT nx, PLINT ny, 00417 PLINT ( *defined )( PLFLT, PLFLT ), 00418 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00419 PLFLT shade_min, PLFLT shade_max, 00420 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00421 PLINT min_color, PLINT min_width, 00422 PLINT max_color, PLINT max_width, 00423 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00424 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00425 PLPointer pltr_data ) 00426 { 00427 plshade_int( zops->f2eval, zp, 00428 NULL, NULL, 00429 /* plc2eval, (PLPointer) &cgrid,*/ 00430 defined, MISSING_MIN_DEF, MISSING_MAX_DEF, nx, ny, xmin, 00431 xmax, ymin, ymax, shade_min, shade_max, 00432 sh_cmap, sh_color, sh_width, 00433 min_color, min_width, max_color, max_width, 00434 fill, rectangular, pltr, pltr_data ); 00435 } 00436 00437 /*----------------------------------------------------------------------*\ 00438 * plshade_int() 00439 * 00440 * Shade region -- this routine does all the work 00441 * 00442 * This routine is internal so the arguments can and will change. 00443 * To retain some compatibility between versions, you must go through 00444 * some stub routine! 00445 * 00446 * 4/95 00447 * 00448 * new: missing_min, missing_max 00449 * 00450 * if data <= missing_max and data >= missing_min 00451 * then the data will beconsidered to be missing 00452 * this allows 2nd way to set undefined points (good for ftn) 00453 * if missing feature is not used, set missing_max < missing_min 00454 * 00455 * parameters: 00456 * 00457 * f2eval, f2eval_data: data to plot 00458 * c2eval, c2eval_data: defined mask (not implimented) 00459 * defined: defined mask (old API - implimented) 00460 * missing_min, missing_max: yet another way to set data to undefined 00461 * nx, ny: array dimensions 00462 * xmin, xmax, ymin, ymax: grid coordinates 00463 * shade_min, shade_max: shade region with values between ... 00464 * sh_cmap, sh_color, sh_width: shading parameters, width is only for hatching 00465 * min_color, min_width: line parameters for boundary (minimum) 00466 * max_color, max_width: line parameters for boundary (maximum) 00467 * set min_width == 0 and max_width == 0 for no contours 00468 * fill: fill function, set to NULL for no shading (contour plot) 00469 * rectangular: flag set to 1 if pltr() maps rectangles to rectangles 00470 * this helps optimize the plotting 00471 * pltr: function to map from grid to plot coordinates 00472 * 00473 * 00474 \*----------------------------------------------------------------------*/ 00475 00476 static void 00477 plshade_int( PLFLT ( *f2eval )( PLINT, PLINT, PLPointer ), 00478 PLPointer f2eval_data, 00479 PLFLT ( *c2eval )( PLINT, PLINT, PLPointer ), 00480 PLPointer c2eval_data, 00481 PLINT ( *defined )( PLFLT, PLFLT ), 00482 PLFLT missing_min, PLFLT missing_max, 00483 PLINT nx, PLINT ny, 00484 PLFLT xmin, PLFLT xmax, PLFLT ymin, PLFLT ymax, 00485 PLFLT shade_min, PLFLT shade_max, 00486 PLINT sh_cmap, PLFLT sh_color, PLINT sh_width, 00487 PLINT min_color, PLINT min_width, 00488 PLINT max_color, PLINT max_width, 00489 void ( *fill )( PLINT, PLFLT *, PLFLT * ), PLINT rectangular, 00490 void ( *pltr )( PLFLT, PLFLT, PLFLT *, PLFLT *, PLPointer ), 00491 PLPointer pltr_data ) 00492 { 00493 PLINT init_width, n, slope = 0, ix, iy; 00494 int count, i, j, nxny; 00495 PLFLT *a, *a0, *a1, dx, dy; 00496 PLFLT x[8], y[8], xp[2], tx, ty; 00497 int *c, *c0, *c1; 00498 00499 if ( plsc->level < 3 ) 00500 { 00501 plabort( "plfshade: window must be set up first" ); 00502 return; 00503 } 00504 00505 if ( nx <= 0 || ny <= 0 ) 00506 { 00507 plabort( "plfshade: nx and ny must be positive" ); 00508 return; 00509 } 00510 00511 if ( shade_min >= shade_max ) 00512 { 00513 plabort( "plfshade: shade_max must exceed shade_min" ); 00514 return; 00515 } 00516 00517 if ( pltr == NULL && plsc->coordinate_transform == NULL ) 00518 rectangular = 1; 00519 00520 int_val = shade_max - shade_min; 00521 init_width = plsc->width; 00522 00523 pen_col_min = min_color; 00524 pen_col_max = max_color; 00525 00526 pen_wd_min = min_width; 00527 pen_wd_max = max_width; 00528 00529 plstyl( (PLINT) 0, NULL, NULL ); 00530 plwid( sh_width ); 00531 if ( fill != NULL ) 00532 { 00533 switch ( sh_cmap ) 00534 { 00535 case 0: 00536 plcol0( (PLINT) sh_color ); 00537 break; 00538 case 1: 00539 plcol1( sh_color ); 00540 break; 00541 default: 00542 plabort( "plfshade: invalid color map selection" ); 00543 return; 00544 } 00545 } 00546 /* alloc space for value array, and initialize */ 00547 /* This is only a temporary kludge */ 00548 nxny = nx * ny; 00549 if ( ( a = (PLFLT *) malloc( nxny * sizeof ( PLFLT ) ) ) == NULL ) 00550 { 00551 plabort( "plfshade: unable to allocate memory for value array" ); 00552 return; 00553 } 00554 00555 for ( ix = 0; ix < nx; ix++ ) 00556 for ( iy = 0; iy < ny; iy++ ) 00557 a[iy + ix * ny] = f2eval( ix, iy, f2eval_data ); 00558 00559 /* alloc space for condition codes */ 00560 00561 if ( ( c = (int *) malloc( nxny * sizeof ( int ) ) ) == NULL ) 00562 { 00563 plabort( "plfshade: unable to allocate memory for condition codes" ); 00564 free( a ); 00565 return; 00566 } 00567 00568 sh_min = shade_min; 00569 sh_max = shade_max; 00570 00571 set_cond( c, a, nxny ); 00572 dx = ( xmax - xmin ) / ( nx - 1 ); 00573 dy = ( ymax - ymin ) / ( ny - 1 ); 00574 a0 = a; 00575 a1 = a + ny; 00576 c0 = c; 00577 c1 = c + ny; 00578 00579 for ( ix = 0; ix < nx - 1; ix++ ) 00580 { 00581 for ( iy = 0; iy < ny - 1; iy++ ) 00582 { 00583 count = c0[iy] + c0[iy + 1] + c1[iy] + c1[iy + 1]; 00584 00585 /* No filling needs to be done for these cases */ 00586 00587 if ( count >= UNDEF ) 00588 continue; 00589 if ( count == 4 * POS ) 00590 continue; 00591 if ( count == 4 * NEG ) 00592 continue; 00593 00594 /* Entire rectangle can be filled */ 00595 00596 if ( count == 4 * OK ) 00597 { 00598 /* find biggest rectangle that fits */ 00599 if ( rectangular ) 00600 { 00601 big_recl( c0 + iy, ny, nx - ix, ny - iy, &i, &j ); 00602 } 00603 else 00604 { 00605 i = j = 1; 00606 } 00607 x[0] = x[1] = ix; 00608 x[2] = x[3] = ix + i; 00609 y[0] = y[3] = iy; 00610 y[1] = y[2] = iy + j; 00611 00612 if ( pltr ) 00613 { 00614 for ( i = 0; i < 4; i++ ) 00615 { 00616 ( *pltr )( x[i], y[i], &tx, &ty, pltr_data ); 00617 x[i] = tx; 00618 y[i] = ty; 00619 } 00620 } 00621 else 00622 { 00623 for ( i = 0; i < 4; i++ ) 00624 { 00625 x[i] = xmin + x[i] * dx; 00626 y[i] = ymin + y[i] * dy; 00627 } 00628 } 00629 if ( fill != NULL ) 00630 exfill( fill, defined, (PLINT) 4, x, y ); 00631 iy += j - 1; 00632 continue; 00633 } 00634 00635 /* Only part of rectangle can be filled */ 00636 00637 n_point = min_points = max_points = 0; 00638 n = find_interval( a0[iy], a0[iy + 1], c0[iy], c0[iy + 1], xp ); 00639 for ( j = 0; j < n; j++ ) 00640 { 00641 x[j] = ix; 00642 y[j] = iy + xp[j]; 00643 } 00644 00645 i = find_interval( a0[iy + 1], a1[iy + 1], 00646 c0[iy + 1], c1[iy + 1], xp ); 00647 00648 for ( j = 0; j < i; j++ ) 00649 { 00650 x[j + n] = ix + xp[j]; 00651 y[j + n] = iy + 1; 00652 } 00653 n += i; 00654 00655 i = find_interval( a1[iy + 1], a1[iy], c1[iy + 1], c1[iy], xp ); 00656 for ( j = 0; j < i; j++ ) 00657 { 00658 x[n + j] = ix + 1; 00659 y[n + j] = iy + 1 - xp[j]; 00660 } 00661 n += i; 00662 00663 i = find_interval( a1[iy], a0[iy], c1[iy], c0[iy], xp ); 00664 for ( j = 0; j < i; j++ ) 00665 { 00666 x[n + j] = ix + 1 - xp[j]; 00667 y[n + j] = iy; 00668 } 00669 n += i; 00670 00671 if ( pltr ) 00672 { 00673 for ( i = 0; i < n; i++ ) 00674 { 00675 ( *pltr )( x[i], y[i], &tx, &ty, pltr_data ); 00676 x[i] = tx; 00677 y[i] = ty; 00678 } 00679 } 00680 else 00681 { 00682 for ( i = 0; i < n; i++ ) 00683 { 00684 x[i] = xmin + x[i] * dx; 00685 y[i] = ymin + y[i] * dy; 00686 } 00687 } 00688 00689 if ( min_points == 4 ) 00690 slope = plctestez( a, nx, ny, ix, iy, shade_min ); 00691 if ( max_points == 4 ) 00692 slope = plctestez( a, nx, ny, ix, iy, shade_max ); 00693 00694 /* n = number of end of line segments */ 00695 /* min_points = number times shade_min meets edge */ 00696 /* max_points = number times shade_max meets edge */ 00697 00698 /* special cases: check number of times a contour is in a box */ 00699 00700 switch ( ( min_points << 3 ) + max_points ) 00701 { 00702 case 000: 00703 case 020: 00704 case 002: 00705 case 022: 00706 if ( fill != NULL && n > 0 ) 00707 exfill( fill, defined, n, x, y ); 00708 break; 00709 case 040: /* 2 contour lines in box */ 00710 case 004: 00711 if ( n != 6 ) 00712 fprintf( stderr, "plfshade err n=%d !6", (int) n ); 00713 if ( slope == 1 && c0[iy] == OK ) 00714 { 00715 if ( fill != NULL ) 00716 exfill( fill, defined, n, x, y ); 00717 } 00718 else if ( slope == 1 ) 00719 { 00720 selected_polygon( fill, defined, x, y, 0, 1, 2, -1 ); 00721 selected_polygon( fill, defined, x, y, 3, 4, 5, -1 ); 00722 } 00723 else if ( c0[iy + 1] == OK ) 00724 { 00725 if ( fill != NULL ) 00726 exfill( fill, defined, n, x, y ); 00727 } 00728 else 00729 { 00730 selected_polygon( fill, defined, x, y, 0, 1, 5, -1 ); 00731 selected_polygon( fill, defined, x, y, 2, 3, 4, -1 ); 00732 } 00733 break; 00734 case 044: 00735 if ( n != 8 ) 00736 fprintf( stderr, "plfshade err n=%d !8", (int) n ); 00737 if ( slope == 1 ) 00738 { 00739 selected_polygon( fill, defined, x, y, 0, 1, 2, 3 ); 00740 selected_polygon( fill, defined, x, y, 4, 5, 6, 7 ); 00741 } 00742 else 00743 { 00744 selected_polygon( fill, defined, x, y, 0, 1, 6, 7 ); 00745 selected_polygon( fill, defined, x, y, 2, 3, 4, 5 ); 00746 } 00747 break; 00748 case 024: 00749 case 042: 00750 /* 3 contours */ 00751 if ( n != 7 ) 00752 fprintf( stderr, "plfshade err n=%d !7", (int) n ); 00753 00754 if ( ( c0[iy] == OK || c1[iy + 1] == OK ) && slope == 1 ) 00755 { 00756 if ( fill != NULL ) 00757 exfill( fill, defined, n, x, y ); 00758 } 00759 else if ( ( c0[iy + 1] == OK || c1[iy] == OK ) && slope == 0 ) 00760 { 00761 if ( fill != NULL ) 00762 exfill( fill, defined, n, x, y ); 00763 } 00764 00765 else if ( c0[iy] == OK ) 00766 { 00767 selected_polygon( fill, defined, x, y, 0, 1, 6, -1 ); 00768 selected_polygon( fill, defined, x, y, 2, 3, 4, 5 ); 00769 } 00770 else if ( c0[iy + 1] == OK ) 00771 { 00772 selected_polygon( fill, defined, x, y, 0, 1, 2, -1 ); 00773 selected_polygon( fill, defined, x, y, 3, 4, 5, 6 ); 00774 } 00775 else if ( c1[iy + 1] == OK ) 00776 { 00777 selected_polygon( fill, defined, x, y, 0, 1, 5, 6 ); 00778 selected_polygon( fill, defined, x, y, 2, 3, 4, -1 ); 00779 } 00780 else if ( c1[iy] == OK ) 00781 { 00782 selected_polygon( fill, defined, x, y, 0, 1, 2, 3 ); 00783 selected_polygon( fill, defined, x, y, 4, 5, 6, -1 ); 00784 } 00785 else 00786 { 00787 fprintf( stderr, "plfshade err logic case 024:042\n" ); 00788 } 00789 break; 00790 default: 00791 fprintf( stderr, "prog err switch\n" ); 00792 break; 00793 } 00794 draw_boundary( slope, x, y ); 00795 00796 if ( fill != NULL ) 00797 { 00798 plwid( sh_width ); 00799 if ( sh_cmap == 0 ) 00800 plcol0( (PLINT) sh_color ); 00801 else if ( sh_cmap == 1 ) 00802 plcol1( sh_color ); 00803 } 00804 } 00805 00806 a0 = a1; 00807 c0 = c1; 00808 a1 += ny; 00809 c1 += ny; 00810 } 00811 00812 free( c ); 00813 free( a ); 00814 plwid( init_width ); 00815 } 00816 00817 /*----------------------------------------------------------------------*\ 00818 * set_cond() 00819 * 00820 * Fills out condition code array. 00821 \*----------------------------------------------------------------------*/ 00822 00823 static void 00824 set_cond( register int *cond, register PLFLT *a, register PLINT n ) 00825 { 00826 while ( n-- ) 00827 { 00828 if ( *a < sh_min ) 00829 *cond++ = NEG; 00830 else if ( *a > sh_max ) 00831 *cond++ = POS; 00832 else 00833 *cond++ = OK; 00834 a++; 00835 } 00836 } 00837 00838 /*----------------------------------------------------------------------*\ 00839 * find_interval() 00840 * 00841 * Two points x(0) = a0, (condition code c0) x(1) = a1, (condition code c1) 00842 * return interval on the line that are shaded 00843 * 00844 * returns 0 : no points to be shaded 1 : x[0] <= x < 1 is the interval 2 : 00845 * x[0] <= x <= x[1] < 1 interval to be shaded n_point, max_points, 00846 * min_points are incremented location of min/max_points are stored 00847 \*----------------------------------------------------------------------*/ 00848 00849 static int 00850 find_interval( PLFLT a0, PLFLT a1, PLINT c0, PLINT c1, PLFLT *x ) 00851 { 00852 register int n; 00853 00854 n = 0; 00855 if ( c0 == OK ) 00856 { 00857 x[n++] = 0.0; 00858 n_point++; 00859 } 00860 if ( c0 == c1 ) 00861 return n; 00862 00863 if ( c0 == NEG || c1 == POS ) 00864 { 00865 if ( c0 == NEG ) 00866 { 00867 x[n++] = linear( a0, a1, sh_min ); 00868 min_pts[min_points++] = n_point++; 00869 } 00870 if ( c1 == POS ) 00871 { 00872 x[n++] = linear( a0, a1, sh_max ); 00873 max_pts[max_points++] = n_point++; 00874 } 00875 } 00876 if ( c0 == POS || c1 == NEG ) 00877 { 00878 if ( c0 == POS ) 00879 { 00880 x[n++] = linear( a0, a1, sh_max ); 00881 max_pts[max_points++] = n_point++; 00882 } 00883 if ( c1 == NEG ) 00884 { 00885 x[n++] = linear( a0, a1, sh_min ); 00886 min_pts[min_points++] = n_point++; 00887 } 00888 } 00889 return n; 00890 } 00891 00892 /*----------------------------------------------------------------------*\ 00893 * selected_polygon() 00894 * 00895 * Draws a polygon from points in x[] and y[]. 00896 * Point selected by v1..v4 00897 \*----------------------------------------------------------------------*/ 00898 00899 static void 00900 selected_polygon( void ( *fill )( PLINT, PLFLT *, PLFLT * ), 00901 PLINT ( *defined )( PLFLT, PLFLT ), 00902 PLFLT *x, PLFLT *y, PLINT v1, PLINT v2, PLINT v3, PLINT v4 ) 00903 { 00904 register PLINT n = 0; 00905 PLFLT xx[4], yy[4]; 00906 00907 if ( fill == NULL ) 00908 return; 00909 if ( v1 >= 0 ) 00910 { 00911 xx[n] = x[v1]; 00912 yy[n++] = y[v1]; 00913 } 00914 if ( v2 >= 0 ) 00915 { 00916 xx[n] = x[v2]; 00917 yy[n++] = y[v2]; 00918 } 00919 if ( v3 >= 0 ) 00920 { 00921 xx[n] = x[v3]; 00922 yy[n++] = y[v3]; 00923 } 00924 if ( v4 >= 0 ) 00925 { 00926 xx[n] = x[v4]; 00927 yy[n++] = y[v4]; 00928 } 00929 exfill( fill, defined, n, xx, yy ); 00930 } 00931 00932 /*----------------------------------------------------------------------*\ 00933 * bisect() 00934 * 00935 * Find boundary recursively by bisection. 00936 * (x1, y1) is in the defined region, while (x2, y2) in the undefined one. 00937 * The result is returned in 00938 \*----------------------------------------------------------------------*/ 00939 00940 static void 00941 bisect( PLINT ( *defined )( PLFLT, PLFLT ), PLINT niter, 00942 PLFLT x1, PLFLT y1, PLFLT x2, PLFLT y2, PLFLT* xb, PLFLT* yb ) 00943 { 00944 PLFLT xm; 00945 PLFLT ym; 00946 00947 if ( niter == 0 ) 00948 { 00949 *xb = x1; 00950 *yb = y1; 00951 return; 00952 } 00953 00954 xm = ( x1 + x2 ) / 2.; 00955 ym = ( y1 + y2 ) / 2.; 00956 00957 if ( defined( xm, ym ) ) 00958 bisect( defined, niter - 1, xm, ym, x2, y2, xb, yb ); 00959 else 00960 bisect( defined, niter - 1, x1, y1, xm, ym, xb, yb ); 00961 } 00962 00963 /*----------------------------------------------------------------------*\ 00964 * exfill() 00965 * 00966 * Fills a polygon from points in x[] and y[] with all points in 00967 * undefined regions dropped and replaced by points at the bisected 00968 * edge of the defined region. 00969 * Note, undefined regions that are confined to the areas between polygon 00970 * points are completely ignored. Also, a range of undefined polygon points 00971 * are simply replaced with a straight line with accurately bisected end 00972 * points. So this routine can produce problematic plotted results 00973 * if the polygon is not a lot smaller than the typical resolution of 00974 * the defined region. 00975 \*----------------------------------------------------------------------*/ 00976 00977 static void 00978 exfill( void ( *fill )( PLINT, PLFLT *, PLFLT * ), 00979 PLINT ( *defined )( PLFLT, PLFLT ), 00980 int n, PLFLT *x, PLFLT *y ) 00981 { 00982 if ( n < 3 ) 00983 { 00984 plabort( "exfill: Not enough points in object" ); 00985 return; 00986 } 00987 00988 if ( defined == NULL ) 00989 00990 ( *fill )( n, x, y ); 00991 00992 else 00993 { 00994 PLFLT *xx; 00995 PLFLT *yy; 00996 PLFLT xb, yb; 00997 PLINT count = 0; 00998 PLINT im1 = n - 1; 00999 PLINT is_defined = defined( x[im1], y[im1] ); 01000 PLINT i; 01001 01002 /* Slightly less than 2 n points are required for xx, yy, but 01003 * allocate room for 2 n to be safe. */ 01004 if ( ( xx = (PLFLT *) malloc( 2 * n * sizeof ( PLFLT ) ) ) == NULL ) 01005 plexit( "exfill: out of memory for xx" ); 01006 if ( ( yy = (PLFLT *) malloc( 2 * n * sizeof ( PLFLT ) ) ) == NULL ) 01007 plexit( "exfill: out of memory for yy." ); 01008 01009 for ( i = 0; i < n; i++ ) 01010 { 01011 /* is_defined tells whether im1 point was in defined region. */ 01012 if ( defined( x[i], y[i] ) ) 01013 { 01014 if ( !is_defined ) 01015 { 01016 /* Cross from undefined (at im1) to defined region. 01017 * Bisect for the first point inside the defined region 01018 * and add it to xx, yy. */ 01019 bisect( defined, NUMBER_BISECTIONS, 01020 x[i], y[i], x[im1], y[im1], &xb, &yb ); 01021 xx[count] = xb; 01022 yy[count++] = yb; 01023 } 01024 /* x[i], y[i] known to be in defined region so add this 01025 * point to xx, yy. */ 01026 xx[count] = x[i]; 01027 yy[count++] = y[i]; 01028 is_defined = 1; 01029 } 01030 else 01031 { 01032 if ( is_defined ) 01033 { 01034 /* Cross from defined (at im1) to undefined region. 01035 * Bisect for the last point in the defined region and 01036 * add it to xx, yy. */ 01037 bisect( defined, NUMBER_BISECTIONS, 01038 x[im1], y[im1], x[i], y[i], &xb, &yb ); 01039 xx[count] = xb; 01040 yy[count++] = yb; 01041 is_defined = 0; 01042 } 01043 } 01044 im1 = i; 01045 } 01046 01047 if ( count >= 3 ) 01048 ( *fill )( count, xx, yy ); 01049 01050 free( xx ); 01051 free( yy ); 01052 } 01053 } 01054 01055 /*----------------------------------------------------------------------*\ 01056 * big_recl() 01057 * 01058 * find a big rectangle for shading 01059 * 01060 * 2 goals: minimize calls to (*fill)() 01061 * keep ratio 1:3 for biggest rectangle 01062 * 01063 * only called by plshade() 01064 * 01065 * assumed that a 1 x 1 square already fits 01066 * 01067 * c[] = condition codes 01068 * ny = c[1][0] == c[ny] (you know what I mean) 01069 * 01070 * returns ix, iy = length of rectangle in grid units 01071 * 01072 * ix < dx - 1 01073 * iy < dy - 1 01074 * 01075 * If iy == 1 -> ix = 1 (so that cond code can be set to skip) 01076 \*----------------------------------------------------------------------*/ 01077 01078 #define RATIO 3 01079 #define COND( x, y ) cond_code[x * ny + y] 01080 01081 static void 01082 big_recl( int *cond_code, register int ny, int dx, int dy, 01083 int *ix, int *iy ) 01084 { 01085 int ok_x, ok_y, j; 01086 register int i, x, y; 01087 register int *cond; 01088 01089 /* ok_x = ok to expand in x direction */ 01090 /* x = current number of points in x direction */ 01091 01092 ok_x = ok_y = 1; 01093 x = y = 2; 01094 01095 while ( ok_x || ok_y ) 01096 { 01097 #ifdef RATIO 01098 if ( RATIO * x <= y || RATIO * y <= x ) 01099 break; 01100 #endif 01101 if ( ok_y ) 01102 { 01103 /* expand in vertical */ 01104 ok_y = 0; 01105 if ( y == dy ) 01106 continue; 01107 cond = &COND( 0, y ); 01108 for ( i = 0; i < x; i++ ) 01109 { 01110 if ( *cond != OK ) 01111 break; 01112 cond += ny; 01113 } 01114 if ( i == x ) 01115 { 01116 /* row is ok */ 01117 y++; 01118 ok_y = 1; 01119 } 01120 } 01121 if ( ok_x ) 01122 { 01123 if ( y == 2 ) 01124 break; 01125 /* expand in x direction */ 01126 ok_x = 0; 01127 if ( x == dx ) 01128 continue; 01129 cond = &COND( x, 0 ); 01130 for ( i = 0; i < y; i++ ) 01131 { 01132 if ( *cond++ != OK ) 01133 break; 01134 } 01135 if ( i == y ) 01136 { 01137 /* column is OK */ 01138 x++; 01139 ok_x = 1; 01140 } 01141 } 01142 } 01143 01144 /* found the largest rectangle of 'ix' by 'iy' */ 01145 *ix = --x; 01146 *iy = --y; 01147 01148 /* set condition code to UNDEF in interior of rectangle */ 01149 01150 for ( i = 1; i < x; i++ ) 01151 { 01152 cond = &COND( i, 1 ); 01153 for ( j = 1; j < y; j++ ) 01154 { 01155 *cond++ = UNDEF; 01156 } 01157 } 01158 } 01159 01160 /*----------------------------------------------------------------------*\ 01161 * draw_boundary() 01162 * 01163 * Draw boundaries of contour regions based on min_pts[], and max_pts[]. 01164 \*----------------------------------------------------------------------*/ 01165 01166 static void 01167 draw_boundary( PLINT slope, PLFLT *x, PLFLT *y ) 01168 { 01169 int i; 01170 01171 if ( pen_col_min != 0 && pen_wd_min != 0 && min_points != 0 ) 01172 { 01173 plcol0( pen_col_min ); 01174 plwid( pen_wd_min ); 01175 if ( min_points == 4 && slope == 0 ) 01176 { 01177 /* swap points 1 and 3 */ 01178 i = min_pts[1]; 01179 min_pts[1] = min_pts[3]; 01180 min_pts[3] = i; 01181 } 01182 pljoin( x[min_pts[0]], y[min_pts[0]], x[min_pts[1]], y[min_pts[1]] ); 01183 if ( min_points == 4 ) 01184 { 01185 pljoin( x[min_pts[2]], y[min_pts[2]], x[min_pts[3]], 01186 y[min_pts[3]] ); 01187 } 01188 } 01189 if ( pen_col_max != 0 && pen_wd_max != 0 && max_points != 0 ) 01190 { 01191 plcol0( pen_col_max ); 01192 plwid( pen_wd_max ); 01193 if ( max_points == 4 && slope == 0 ) 01194 { 01195 /* swap points 1 and 3 */ 01196 i = max_pts[1]; 01197 max_pts[1] = max_pts[3]; 01198 max_pts[3] = i; 01199 } 01200 pljoin( x[max_pts[0]], y[max_pts[0]], x[max_pts[1]], y[max_pts[1]] ); 01201 if ( max_points == 4 ) 01202 { 01203 pljoin( x[max_pts[2]], y[max_pts[2]], x[max_pts[3]], 01204 y[max_pts[3]] ); 01205 } 01206 } 01207 } 01208 01209 /*----------------------------------------------------------------------*\ 01210 * 01211 * plctest( &(x[0][0]), PLFLT level) 01212 * where x was defined as PLFLT x[4][4]; 01213 * 01214 * determines if the contours associated with level have 01215 * positive slope or negative slope in the box: 01216 * 01217 * (2,3) (3,3) 01218 * 01219 * (2,2) (3,2) 01220 * 01221 * this is heuristic and can be changed by the user 01222 * 01223 * return 1 if positive slope 01224 * 0 if negative slope 01225 * 01226 * algorithmn: 01227 * 1st test: 01228 * find length of contours assuming positive and negative slopes 01229 * if the length of the negative slope contours is much bigger 01230 * than the positive slope, then the slope is positive. 01231 * (and vice versa) 01232 * (this test tries to minimize the length of contours) 01233 * 01234 * 2nd test: 01235 * if abs((top-right corner) - (bottom left corner)) > 01236 * abs((top-left corner) - (bottom right corner)) ) then 01237 * return negatiave slope. 01238 * (this test tries to keep the slope for different contour levels 01239 * the same) 01240 \*----------------------------------------------------------------------*/ 01241 01242 #define X( a, b ) ( x[a * 4 + b] ) 01243 #define POSITIVE_SLOPE (PLINT) 1 01244 #define NEGATIVE_SLOPE (PLINT) 0 01245 #define RATIO_SQ 6.0 01246 01247 static PLINT 01248 plctest( PLFLT *x, PLFLT level ) 01249 { 01250 int i, j; 01251 double t[4], sorted[4], temp; 01252 01253 sorted[0] = t[0] = X( 1, 1 ); 01254 sorted[1] = t[1] = X( 2, 2 ); 01255 sorted[2] = t[2] = X( 1, 2 ); 01256 sorted[3] = t[3] = X( 2, 1 ); 01257 01258 for ( j = 1; j < 4; j++ ) 01259 { 01260 temp = sorted[j]; 01261 i = j - 1; 01262 while ( i >= 0 && sorted[i] > temp ) 01263 { 01264 sorted[i + 1] = sorted[i]; 01265 i--; 01266 } 01267 sorted[i + 1] = temp; 01268 } 01269 /* sorted[0] == min */ 01270 01271 /* find min contour */ 01272 temp = int_val * ceil( sorted[0] / int_val ); 01273 if ( temp < sorted[1] ) 01274 { 01275 /* one contour line */ 01276 for ( i = 0; i < 4; i++ ) 01277 { 01278 if ( t[i] < temp ) 01279 return i / 2; 01280 } 01281 } 01282 01283 /* find max contour */ 01284 temp = int_val * floor( sorted[3] / int_val ); 01285 if ( temp > sorted[2] ) 01286 { 01287 /* one contour line */ 01288 for ( i = 0; i < 4; i++ ) 01289 { 01290 if ( t[i] > temp ) 01291 return i / 2; 01292 } 01293 } 01294 /* nothing better to do - be consistant */ 01295 return POSITIVE_SLOPE; 01296 } 01297 01298 /*----------------------------------------------------------------------*\ 01299 * plctestez 01300 * 01301 * second routine - easier to use 01302 * fills in x[4][4] and calls plctest 01303 * 01304 * test location a[ix][iy] (lower left corner) 01305 \*----------------------------------------------------------------------*/ 01306 01307 static PLINT 01308 plctestez( PLFLT *a, PLINT nx, PLINT ny, PLINT ix, 01309 PLINT iy, PLFLT level ) 01310 { 01311 PLFLT x[4][4]; 01312 int i, j, ii, jj; 01313 01314 for ( i = 0; i < 4; i++ ) 01315 { 01316 ii = ix + i - 1; 01317 ii = MAX( 0, ii ); 01318 ii = MIN( ii, nx - 1 ); 01319 for ( j = 0; j < 4; j++ ) 01320 { 01321 jj = iy + j - 1; 01322 jj = MAX( 0, jj ); 01323 jj = MIN( jj, ny - 1 ); 01324 x[i][j] = a[ii * ny + jj]; 01325 } 01326 } 01327 return plctest( &( x[0][0] ), level ); 01328 }