/*@@ @file Write2D.c @date Thu May 11 2000 @author Thomas Radke @desc Output two-dimensional slices in ASCII gnuplot format. @enddesc @@*/ #include #include #include #include "cctk.h" #include "cctk_Parameters.h" #include "CactusPUGH/PUGHSlab/src/PUGHSlab.h" #include "CactusBase/IOUtil/src/ioGH.h" #include "ioASCIIGH.h" /* enable debug output */ /*#define IOASCII_DEBUG*/ /* macro to output a time slice as typed data */ #define OUTPUT_TYPED_DATA(coord_hsizes, coord_i, coord_j, stagger_offset_i, \ stagger_offset_j, data_hsizes, cctk_type, c_type, \ data, fmt_string, file) \ { \ int i, j; \ cctk_type *typed_data = (cctk_type *) data; \ \ \ for (j = 0; j < data_hsizes [1]; j++) \ { \ for (i = 0; i < data_hsizes [0]; i++) \ { \ fprintf (file, fmt_string, \ (double) (*coord_i++ + stagger_offset_i), \ (double) (*coord_j++ + stagger_offset_j), \ (c_type) *typed_data++); \ } \ fprintf (file, "\n"); \ coord_i += coord_hsizes [0] - data_hsizes [0]; \ coord_j += coord_hsizes [0] - data_hsizes [0]; \ } \ coord_i -= coord_hsizes [0] * data_hsizes [1]; \ coord_j -= coord_hsizes [0] * data_hsizes [1]; \ } /*@@ @routine IOASCII_Write2D @date Thu May 11 2000 @author Thomas Radke @desc Writes the 2D slices of a variable into separate ASCII files @enddesc @calledby IOASCII_Output2DVarAs IOASCII_TriggerOutput2D @history @endhistory @var GH @vdesc Pointer to CCTK GH @vtype cGH @vio in @vcomment @endvar @var index @vdesc index of variable to output @vtype int @vio in @vcomment @endvar @var alias @vdesc alias name of variable to output @vtype const char * @vio in @vcomment @endvar @@*/ void IOASCII_Write2D (cGH *GH, int index, const char *alias) { DECLARE_CCTK_PARAMETERS int myproc; ioGH *ioUtilGH; asciiioGH *myGH; const char *header_fmt_string; /* header format string */ const char *data_fmt_string_int; /* data format string for int types */ const char *data_fmt_string_real; /* data format string for float types */ int dir; int timelevel; int groupindex; cGroup groupinfo; FILE **fdset_2D; /* array of output file pointers */ int coord_index [3]; /* variable indices for xyz coordinates */ int coord_timelevel [3]; /* coordinates' current timelevels */ int origin [3]; /* the slice origin */ const char *errormsg; /* error message string */ /* to make the compiler happy */ fdset_2D = NULL; /* get the variable group indormation */ groupindex = CCTK_GroupIndexFromVarI (index); CCTK_GroupData (groupindex, &groupinfo); /* check if variable is of dimension 3 and has storage assigned */ errormsg = NULL; if (groupinfo.dim != 3) errormsg = "No IOASCII 2D output for '%s' (dim != 3)"; else if (! CCTK_QueryGroupStorageI (GH, groupindex)) errormsg = "No IOASCII 2D output for '%s' (no storage)"; if (errormsg) { char *fullname = CCTK_FullName (index); CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING, errormsg, fullname); free (fullname); return; } /* Get the handles for IOUtil and IOASCII extensions */ ioUtilGH = (ioGH *) GH->extensions [CCTK_GHExtensionHandle ("IO")]; myGH = (asciiioGH *) GH->extensions [CCTK_GHExtensionHandle ("IOASCII")]; /* set header and data format strings */ if (CCTK_Equals (out_format, "f")) { header_fmt_string = "\n\n#Time = %f\n"; data_fmt_string_int = "%.13f\t\t%.13f\t\t%d\n"; data_fmt_string_real = "%.13f\t\t%.13f\t\t%.13f\n"; } else { header_fmt_string = "\n\n#Time = %e\n"; data_fmt_string_int = "%.13e\t\t%.13e\t\t%d\n"; data_fmt_string_real = "%.13e\t\t%.13e\t\t%.13e\n"; } /* What processor are we on? */ myproc = CCTK_MyProc (GH); /* Open the files on the first trip through if we are proc. 0 */ if (myproc == 0) { /* see if output file for this alias name was already created */ fdset_2D = (FILE **) GetNamedData (myGH->fileList_2D, alias); if (fdset_2D == NULL) { char *fname; fdset_2D = (FILE **) malloc (3 * sizeof (FILE *)); fname = (char *) malloc (strlen (myGH->outdir2D) + strlen (alias) + #if 0 FIXME: get rid of PUGH here strlen (pughGH->identity_string) + 20); #else 20); #endif /* Open/Create files */ for (dir = 0; dir < 3; dir++) { IOUtil_AdvertisedFileDesc_t advertised_file; const char *extensions [3] = {"yz", "xz", "xy"}; #if 0 sprintf (fname, "%s/%s_2d_%s%s.gnuplot", myGH->outdir2D, alias, extensions [dir], pughGH->identity_string); #else sprintf (fname, "%s/%s_2d_%s.gnuplot", myGH->outdir2D, alias, extensions [dir]); #endif /* if restart from recovery, try to open an existing file ... */ fdset_2D [dir] = NULL; if (ioUtilGH->recovered) fdset_2D [dir] = fopen (fname, "a"); /* otherwise or if that failed, create a new one */ if (! fdset_2D [dir]) fdset_2D [dir] = fopen (fname, "w"); if (! fdset_2D [dir]) { CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Cannot open 2D output file '%s'", fname); return; } /* advertise the file for downloading */ advertised_file.slice = (char *) extensions [dir]; advertised_file.thorn = CCTK_THORNSTRING; advertised_file.varname = CCTK_FullName (index); advertised_file.description = "Two-dimensional slice plots"; advertised_file.mimetype = "application/gnuplot"; IOUtil_AdvertiseFile (GH, fname, &advertised_file); free (advertised_file.varname); } /* store file desriptors in database */ StoreNamedData (&myGH->fileList_2D, alias, fdset_2D); free (fname); } } /* get the coordinate indices */ switch (groupinfo.dim) { case 1: coord_index [0] = CCTK_CoordIndex (1,NULL,"cart1d"); break; case 2: coord_index [0] = CCTK_CoordIndex (1,NULL,"cart2d"); coord_index [1] = CCTK_CoordIndex (2,NULL,"cart2d"); case 3: coord_index [0] = CCTK_CoordIndex (1,NULL,"cart3d"); coord_index [1] = CCTK_CoordIndex (2,NULL,"cart3d"); coord_index [2] = CCTK_CoordIndex (3,NULL,"cart3d"); break; default: CCTK_VWarn(4,__LINE__,__FILE__,"IOASCII", "Cannot find appropriate coordinate system"); break; } /* get the coordinate timelevels */ for (dir = 0; dir < 3; dir++) { coord_timelevel [dir] = CCTK_NumTimeLevelsFromVarI (coord_index [dir]) - 1; if (coord_timelevel [dir] > 0) coord_timelevel [dir]--; } /* get the timelevel for the variable to output */ timelevel = CCTK_NumTimeLevelsFromVarI (index) - 1; if (timelevel > 0) timelevel--; for (dir = 0; dir < 3; dir++) { int dir_i, dir_j; int directions [3]; const int lengths [2] = {-1, -1}; const int downsamples [2] = {1, 1}; int coord_hsizes [2], data_hsizes [2]; CCTK_REAL *coord_data_i, *coord_data_j; void *data; /* get the directions to span the hyperslab */ if (dir == 0) { dir_i = 1; dir_j = 2; /* yz */ } else if (dir == 1) { dir_i = 0; dir_j = 2; /* xz */ } else { dir_i = 0; dir_j = 1; /* xy */ } /* set the origin using the slice center from IOUtil */ memset (origin, 0, GH->cctk_dim * sizeof (int)); origin [dir] = myGH->sp2xyz [groupinfo.dim-1][dir]; /* set the directions vector */ memset (directions, 0, sizeof (directions)); directions [dir] = 1; /* get the i-coordinate slice */ if (Hyperslab_GetHyperslab (GH, 0, coord_index [dir_i], coord_timelevel [dir_i], 2, origin, directions, lengths, downsamples, (void **) &coord_data_i, coord_hsizes) < 0) { char *fullname = CCTK_FullName (coord_index [dir_i]); CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Failed to extract 2D hyperslab for variable '%s'", fullname); free (fullname); return; } /* get the j-coordinate slice */ if (Hyperslab_GetHyperslab (GH, 0, coord_index [dir_j], coord_timelevel [dir_j], 2, origin, directions, lengths, downsamples, (void **) &coord_data_j, coord_hsizes) < 0) { char *fullname = CCTK_FullName (coord_index [dir_j]); CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Failed to extract 2D hyperslab for variable '%s'", fullname); free (fullname); free (coord_data_i); return; } /* get the variable slice */ if (Hyperslab_GetHyperslab (GH, 0, index, timelevel, 2, origin, directions, lengths, downsamples, &data, data_hsizes) < 0) { char *fullname = CCTK_FullName (index); CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, "Failed to extract 2D hyperslab for variable '%s'", fullname); free (fullname); free (coord_data_i); free (coord_data_j); return; } /* proc 0 writes */ if (myproc == 0) { CCTK_REAL stagger_offset_i, stagger_offset_j; /* get the staggering offset for the coordinates */ stagger_offset_i = CCTK_StaggerDirIndex (dir_i, groupinfo.stagtype) * 0.5 * GH->cctk_delta_space [dir_i]; stagger_offset_j = CCTK_StaggerDirIndex (dir_j, groupinfo.stagtype) * 0.5 * GH->cctk_delta_space [dir_j]; /* print out header */ fprintf (fdset_2D [dir], header_fmt_string, GH->cctk_time); switch (groupinfo.vartype) { case CCTK_VARIABLE_CHAR: OUTPUT_TYPED_DATA (coord_hsizes, coord_data_i, coord_data_j, stagger_offset_i, stagger_offset_j, data_hsizes, CCTK_CHAR, int, data, data_fmt_string_int, fdset_2D [dir]); break; case CCTK_VARIABLE_INT: OUTPUT_TYPED_DATA (coord_hsizes, coord_data_i, coord_data_j, stagger_offset_i, stagger_offset_j, data_hsizes, CCTK_INT, int, data, data_fmt_string_int, fdset_2D [dir]); break; case CCTK_VARIABLE_REAL: OUTPUT_TYPED_DATA (coord_hsizes, coord_data_i, coord_data_j, stagger_offset_i, stagger_offset_j, data_hsizes, CCTK_REAL, double, data, data_fmt_string_real, fdset_2D [dir]); break; default: CCTK_WARN (1, "Unsupported variable type"); break; } /* keep the file open but flush it */ fflush (fdset_2D [dir]); free (data); free (coord_data_j); free (coord_data_i); } /* end of outputting the data by processor 0 */ } /* end of looping through xyz directions */ }