From 14643128af1ed216faab10b8cf3f939d47242715 Mon Sep 17 00:00:00 2001 From: tradke Date: Thu, 4 Nov 2004 19:04:51 +0000 Subject: Added a steerable integer parameter IOJpeg::refinement_factor which lets you enlarge the resulting JPEG images by some factor. The default is 1 (meaning no resizing, images will have the same size as underlying grid). git-svn-id: http://svn.cactuscode.org/arrangements/CactusIO/IOJpeg/trunk@106 eff87b29-5268-4891-90a3-a07138403961 --- doc/documentation.tex | 4 +++ param.ccl | 5 +++ src/Write.c | 94 +++++++++++++++++++++++++++++++++++++++++++++++++-- 3 files changed, 101 insertions(+), 2 deletions(-) diff --git a/doc/documentation.tex b/doc/documentation.tex index eabdf2b..6feeeef 100644 --- a/doc/documentation.tex +++ b/doc/documentation.tex @@ -40,6 +40,10 @@ the images to be advertised to the Web Server the parameter {\tt iojpeg::mode = Note that the parameter {\tt iojpeg::mode} determines whether a jpeg image is ccreated and kept for individual timesteps, or whether only the image from the current data is kept. Only in the second case is the jpeg file {\tt advertised} (otherwise there are potentially thousands of files advertised to e.g. the {\tt HTTPD} thorn). Also note that using the standard mode and creating jpegs every iteration can quite quickly lead to {\it inode} problems on a machine. +The steerable parameter {\tt IOJpeg::refinement\_factor} determines whether the +resulting JPEG images (of same size as the underlying grid) should be refined +by a certain factor. If refinement is enabled ({\tt IOJpeg::refinement\_factor} $>$ 1) an interpolation will be done to enlarge the images to the requested size. For this case, a thorn providing local interpolation operators must be activated (eg. thorn {\tt CactusBase/LocalInterp}). + We are planning to develop this thorn more to provide more features (eg a range of 2D images from a 3D dataset, add possibility to save images for a movie and advertise current image). diff --git a/param.ccl b/param.ccl index c64ef3a..15d8e89 100644 --- a/param.ccl +++ b/param.ccl @@ -71,6 +71,11 @@ REAL colormap_max "maximum value to be mapped to colors" STEERABLE = ALWAYS *:* :: "Only for custom colormap scale" } +1.0 +INT refinement_factor "Refine each 2D slice by a certain factor (using interpolation) ?" STEERABLE = ALWAYS +{ + 1:* :: "A factor greater 0" +} 1 + ################################ # Choosing what planes to output diff --git a/src/Write.c b/src/Write.c index a4e1612..8c416c9 100644 --- a/src/Write.c +++ b/src/Write.c @@ -28,6 +28,7 @@ CCTK_FILEVERSION(CactusIO_IOJpeg_Write_c) static void WriteData (const cGH *GH, int vindex, const char *alias, int dim, int dir, CCTK_REAL min, CCTK_REAL max, const CCTK_INT hsize[2], const void *hdata); +static void *RefineData (CCTK_INT input_dims[2], const void *input_data); /*@@ @@ -77,7 +78,7 @@ int IOJpeg_Write (const cGH *GH, CCTK_INT vindex, const char *alias) char *fullname; int extent_int[3]; CCTK_INT origin[3], extent[2], direction[6], hsize[2]; - void *hdata; + void *hdata, *outdata; CCTK_REAL min, max; const CCTK_INT htype = CCTK_VARIABLE_REAL; DECLARE_CCTK_PARAMETERS @@ -193,7 +194,18 @@ int IOJpeg_Write (const cGH *GH, CCTK_INT vindex, const char *alias) { if (i == 0) { - WriteData (GH, vindex, alias, gdata.dim, dir, min, max, hsize, hdata); + /* refine the 2D slice if requested */ + outdata = refinement_factor > 1 ? RefineData (hsize, hdata) : hdata; + if (! outdata) + { + CCTK_WARN (1, "Failed to refine 2D array"); + outdata = hdata; + } + WriteData (GH, vindex, alias, gdata.dim, dir, min, max, hsize, outdata); + if (outdata != hdata) + { + free (outdata); + } } /* clean up */ @@ -210,6 +222,84 @@ int IOJpeg_Write (const cGH *GH, CCTK_INT vindex, const char *alias) /******************************************************************** ******************** Internal Routines ************************ ********************************************************************/ +/*@@ + @routine RefineData + @date Thu 4 November 2004 + @author Thomas Radke + @desc + Refines a given 2D input array by a certain factor + @enddesc + + @returntype void * + @returndesc + a pointer to the allocated refined output array (must be freed + by the caller), or NULL in case of an error + @endreturndesc + @@*/ +static void *RefineData (CCTK_INT input_dims[2], const void *input_data) +{ + const CCTK_REAL coord_origin[2] = {0.0, 0.0}; + CCTK_REAL coord_delta[2]; + const CCTK_INT array_type_codes[2] = {CCTK_VARIABLE_REAL, CCTK_VARIABLE_REAL}; + CCTK_REAL *interp_coords[2], *refined_data; + int i, j, interp_handle, N_interp_points; + CCTK_INT output_dims[2]; + DECLARE_CCTK_PARAMETERS + + + interp_handle = CCTK_InterpHandle ("uniform cartesian"); + if (interp_handle < 0) + { + CCTK_WARN (1, "Couldn't get handle for interpolation operator 'uniform " + "cartesian'. Did you forget to activate a thorn providing " + "CCTK_InterpLocalUniform() ?"); + return (NULL); + } + + coord_delta[0] = coord_delta[1] = refinement_factor; + + output_dims[0] = refinement_factor * (input_dims[0] - 1); + output_dims[1] = refinement_factor * (input_dims[1] - 1); + N_interp_points = output_dims[0] * output_dims[1]; + + interp_coords[0] = malloc (2 * N_interp_points * sizeof (CCTK_REAL)); + interp_coords[1] = interp_coords[0]; + for (i = 0; i < output_dims[0]; i++) + { + for (j = 0; j < output_dims[1]; j++) + { + *interp_coords[0]++ = i; + *interp_coords[1]++ = j; + } + } + interp_coords[0] -= N_interp_points; + interp_coords[1] -= N_interp_points; + + refined_data = malloc (N_interp_points * sizeof (CCTK_REAL)); + if (CCTK_InterpLocalUniform (2, interp_handle, -1, + coord_origin, coord_delta, N_interp_points, + CCTK_VARIABLE_REAL, + (const void *const *) interp_coords, + 1, input_dims, array_type_codes, + &input_data, 1, array_type_codes, + (void *const *) &refined_data) != 0) + { + CCTK_WARN (1, "Failed to interpolate 2D array"); + free (refined_data); + refined_data = NULL; + } + else + { + input_dims[0] = output_dims[0]; + input_dims[1] = output_dims[1]; + } + + free (interp_coords[0]); + + return (refined_data); +} + + /*@@ @routine WriteData @date Thu 18 April 2002 -- cgit v1.2.3