From 387489fc876a1cb20c60cfb00a28430320cf2ef3 Mon Sep 17 00:00:00 2001 From: tradke Date: Thu, 12 Jun 2003 12:42:38 +0000 Subject: Fixed IOUtil_1DLines() and IOUtil_2DPlanes() to setup slice centers even if no coordinate system is available but slice index coordinates were specified by the user instead. This closes PR CactusBase/1536. git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/IOUtil/trunk@190 b32723a9-ab3a-4a60-88e2-2e5d99d7c17a --- src/Utils.c | 210 +++++++++++++++++++++++++++++++++++------------------------- 1 file changed, 124 insertions(+), 86 deletions(-) diff --git a/src/Utils.c b/src/Utils.c index 4b559f6..4870d61 100644 --- a/src/Utils.c +++ b/src/Utils.c @@ -301,105 +301,122 @@ int IOUtil_1DLines (const cGH *GH, CCTK_REAL *const *const origin_phys, int *const *slice_center) { - int dim, dir; + int dim, dir, handle, no_origin_index; char coord_system_name[20]; CCTK_REAL *lower_range, *upper_range; - int retval; - retval = 0; - - /* allocate arrays for ranges */ - lower_range = (CCTK_REAL *) calloc (2 * num_dims, sizeof (CCTK_REAL)); - upper_range = lower_range + num_dims; /* get the appropriate coordinate system name */ sprintf (coord_system_name, "cart%dd", num_dims); - if (CCTK_CoordSystemHandle (coord_system_name) >= 0) + handle = CCTK_CoordSystemHandle (coord_system_name); + + /* if no coordinate system is available then origin_index[] must be given */ + no_origin_index = origin_index == NULL; + if (handle < 0 && origin_index) { - /* get the ranges in every direction */ + /* check that origin_index[] is valid */ for (dir = 0; dir < num_dims; dir++) { - if (CCTK_CoordRange (GH, &lower_range[dir], &upper_range[dir], - dir + 1, NULL, coord_system_name) < 0) + for (dim = 0; dim < num_dims; dim++) { - CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_1DLines: Could not get ranges for %c-direction " - "of coordinate system '%s'", - 'x' + dir, coord_system_name); + if (origin_index[dir][dim] < 0) + { + no_origin_index = 1; + break; + } } } } - else + + if (handle < 0 && no_origin_index) { - CCTK_VWarn (4, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_1DLines: Cartesian coordinate system '%s' not found" - " - no slice center set up for output of 1D lines from %dD " - "variables", coord_system_name, num_dims); + CCTK_VWarn (3, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_1DLines: Cartesian coordinate system '%s' not found " + "and no slice center index coordinates set - slice center will " + "not be set up for output of 1D lines from %dD variables", + coord_system_name, num_dims); for (dir = 0; dir < num_dims; dir++) { - for (dim = 0; dim < num_dims; dim++) + memset (slice_center[dir], 0, num_dims * sizeof (int)); + } + + return (-1); + } + + lower_range = upper_range = NULL; + if (handle >= 0) + { + /* allocate arrays for ranges */ + lower_range = calloc (2 * num_dims, sizeof (CCTK_REAL)); + upper_range = lower_range + num_dims; + + /* get the ranges in every direction */ + for (dir = 0; dir < num_dims; dir++) + { + if (CCTK_CoordRange (GH, &lower_range[dir], &upper_range[dir], + dir + 1, NULL, coord_system_name) < 0) { - slice_center[dir][dim] = 0; + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_1DLines: Could not get ranges for %c-direction " + "of coordinate system '%s'", + 'x' + dir, coord_system_name); } } - - retval = -1; } /* now set the slice center for each line according to origin_index[] or origin_phys[] */ - - if (! retval) + for (dir = 0; dir < num_dims; dir++) { - for (dir = 0; dir < num_dims; dir++) + for (dim = 0; dim < num_dims; dim++) { - for (dim = 0; dim < num_dims; dim++) + if (dim == dir) { - if (dim == dir) - { - /* line always starts at the first point */ - slice_center[dir][dim] = 0; - } - else if (origin_index && origin_index[dir][dim] >= 0) - { - /* FIXME: check upper index bounds also ?? */ - slice_center[dir][dim] = origin_index[dir][dim]; - } - else if (lower_range[dim] > origin_phys[dir][dim] || - upper_range[dim] < origin_phys[dir][dim]) - { - CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_1DLines: %c-coordinate (%f) for slice center of " - "1D lines in %c-direction for %dD variables is out of " - "grid coordinates range (%f, %f)", - 'x' + dim, (double) origin_phys[dir][dim], - 'x' + dir, num_dims, - (double) lower_range[dim], (double) upper_range[dim]); - CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_1DLines: slice center will default to %c-index 0", - 'x' + dir); - slice_center[dir][dim] = 0; - } - else - { - /* Find index for first point above the chosen coordinate */ - slice_center[dir][dim] = - ceil ((origin_phys[dir][dim] - lower_range[dim]) / - GH->cctk_delta_space[dim] - 1e-6); + /* line always starts at the first point */ + slice_center[dir][dim] = 0; + } + else if (origin_index && origin_index[dir][dim] >= 0) + { + /* FIXME: check upper index bounds also ?? */ + slice_center[dir][dim] = origin_index[dir][dim]; + } + else if (lower_range[dim] > origin_phys[dir][dim] || + upper_range[dim] < origin_phys[dir][dim]) + { + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_1DLines: %c-coordinate (%f) for slice center of " + "1D lines in %c-direction for %dD variables is out of " + "grid coordinates range (%f, %f)", + 'x' + dim, (double) origin_phys[dir][dim], + 'x' + dir, num_dims, + (double) lower_range[dim], (double) upper_range[dim]); + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_1DLines: slice center will default to %c-index 0", + 'x' + dir); + slice_center[dir][dim] = 0; + } + else + { + /* Find index for first point above the chosen coordinate */ + slice_center[dir][dim] = + ceil ((origin_phys[dir][dim] - lower_range[dim]) / + GH->cctk_delta_space[dim] - 1e-6); #ifdef DEBUG_IOUTIL - printf("spxyz for %c-coord of lines in %c-direction is %d\n", - 'x' + dim,'x' + dir, slice_center[dir][dim]); + printf("spxyz for %c-coord of lines in %c-direction is %d\n", + 'x' + dim,'x' + dir, slice_center[dir][dim]); #endif - } } } } /* free allocated resources */ - free (lower_range); + if (lower_range) + { + free (lower_range); + } - return (retval); + return (0); } @@ -452,37 +469,57 @@ int IOUtil_2DPlanes (const cGH *GH, const CCTK_REAL *origin_phys, int *slice_center) { - int dir; + int handle, dir; char coord_system_name[20]; CCTK_REAL *lower_range, *upper_range; /* get the appropriate coordinate system name */ sprintf (coord_system_name, "cart%dd", num_dims); + handle = CCTK_CoordSystemHandle (coord_system_name); - if (CCTK_CoordSystemHandle (coord_system_name) < 0) + /* if no coordinate system is available then origin_index[] must be given */ + if (handle < 0 && origin_index) { - CCTK_VWarn (4, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_2DPlanes: Cartesian coordinate system '%s' not found" - " - no slice center set up for output of 2D planes from %dD " - "variables", coord_system_name, num_dims); - return (-1); + /* check that origin_index[] is valid */ + for (dir = 0; dir < num_dims; dir++) + { + if (origin_index[dir] < 0) + { + origin_index = NULL; + break; + } + } } - /* allocate arrays for ranges */ - lower_range = (CCTK_REAL *) calloc (2 * num_dims, sizeof (CCTK_REAL)); - upper_range = lower_range + num_dims; + if (handle < 0 && ! origin_index) + { + CCTK_VWarn (3, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_2DPlanes: Cartesian coordinate system '%s' not found " + "and no slice center index coordinates set - slice center will " + "not be set up for output of 2D planes from %dD variables", + coord_system_name, num_dims); + return (-1); + } - /* get the ranges in every direction */ - for (dir = 0; dir < num_dims; dir++) + lower_range = upper_range = NULL; + if (handle >= 0) { - if (CCTK_CoordRange (GH, &lower_range[dir], &upper_range[dir], - num_dims - dir, NULL, coord_system_name) < 0) + /* allocate arrays for ranges */ + lower_range = calloc (2 * num_dims, sizeof (CCTK_REAL)); + upper_range = lower_range + num_dims; + + /* get the ranges in every direction */ + for (dir = 0; dir < num_dims; dir++) { - CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, - "IOUtil_2DPlanes: Could not get ranges for %c-direction " - "of coordinate system '%s'", - 'x' + (num_dims - dir - 1), coord_system_name); + if (CCTK_CoordRange (GH, &lower_range[dir], &upper_range[dir], + num_dims - dir, NULL, coord_system_name) < 0) + { + CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, + "IOUtil_2DPlanes: Could not get ranges for %c-direction " + "of coordinate system '%s'", + 'x' + (num_dims - dir - 1), coord_system_name); + } } } @@ -522,7 +559,10 @@ int IOUtil_2DPlanes (const cGH *GH, } /* free allocated resources */ - free (lower_range); + if (lower_range) + { + free (lower_range); + } return (0); } @@ -704,13 +744,11 @@ static void SetOutputVar (int vindex, const char *optstring, void *arg) char key[128]; CCTK_INT type, nelems; char *fullname; - const info_t *info; + const info_t *info = arg; ioRequest *request; DECLARE_CCTK_PARAMETERS - info = (const info_t *) arg; - /* allocate a new I/O request structure and initialize it with defaults */ request = IOUtil_DefaultIORequest (info->GH, vindex, info->out_every_default); -- cgit v1.2.3