1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
|
/*@@
@file ChooseOutput.c
@author Gabrielle Allen
@date July 6 2000
@desc
Choose what 1D slices and 2D planes to output by IOASCII.
@enddesc
@version $Id$
@@*/
#include <stdlib.h>
#include <string.h>
#include "cctk.h"
#include "cctk_Arguments.h"
#include "cctk_Parameters.h"
#include "CactusBase/IOUtil/src/ioutil_Utils.h"
#include "ioASCIIGH.h"
/* the rcs ID and its dummy function to use it */
static const char *rcsid = "$Header$";
CCTK_FILEVERSION(CactusBase_IOASCII_ChooseOutput_c)
/********************************************************************
******************** Macro Definitions ************************
********************************************************************/
/* macro to choose origin according actual parameter settings:
1. Indices from IOASCII
2. Indices from IOUtil
3. Coords from IOASCII
4. Coords from IOUtil
*/
#define GET_SLICE(IOASCII_param, IOUtil_param, index, coord) \
{ \
index = IOASCII_param##i >= 0 ? IOASCII_param##i : IOUtil_param##i; \
coord = IOASCII_param != -424242 ? IOASCII_param : IOUtil_param; \
}
/********************************************************************
******************** External Routines ************************
********************************************************************/
void IOASCII_Choose1D (CCTK_ARGUMENTS);
void IOASCII_Choose2D (CCTK_ARGUMENTS);
/*@@
@routine IOASCII_Choose1D
@author Gabrielle Allen
@date July 6 2000
@desc
Use parameters to choose the 1D slices through the output data.
@enddesc
@calls IOUtil_1DLines
@var GH
@vdesc pointer to CCTK grid hierarchy
@vtype const cGH *
@vio in
@endvar
@@*/
void IOASCII_Choose1D (CCTK_ARGUMENTS)
{
int i, j, maxdim;
asciiioGH *myGH;
int *origin_index[3];
CCTK_REAL *origin_phys[3];
DECLARE_CCTK_PARAMETERS
/* allocate arrays for origins */
origin_phys[0] = calloc (3 * 3, sizeof (CCTK_REAL));
origin_phys[1] = origin_phys[0] + 3;
origin_phys[2] = origin_phys[1] + 3;
origin_index[0] = calloc (3 * 3, sizeof (int));
origin_index[1] = origin_index[0] + 3;
origin_index[2] = origin_index[1] + 3;
/* get slice points */
GET_SLICE (out1D_xline_y, out_xline_y, origin_index[0][1], origin_phys[0][1]);
GET_SLICE (out1D_xline_z, out_xline_z, origin_index[0][2], origin_phys[0][2]);
GET_SLICE (out1D_yline_x, out_yline_x, origin_index[1][0], origin_phys[1][0]);
GET_SLICE (out1D_yline_z, out_yline_z, origin_index[1][2], origin_phys[1][2]);
GET_SLICE (out1D_zline_x, out_zline_x, origin_index[2][0], origin_phys[2][0]);
GET_SLICE (out1D_zline_y, out_zline_y, origin_index[2][1], origin_phys[2][1]);
maxdim = CCTK_MaxDim ();
myGH = CCTK_GHExtension (cctkGH, "IOASCII");
myGH->spxyz = malloc (maxdim * sizeof (int **));
for (i = 0; i < maxdim; i++)
{
myGH->spxyz[i] = malloc ((i + 1) * sizeof (int *));
for (j = 0; j <= i; j++)
{
myGH->spxyz[i][j] = calloc (i + 1, sizeof (int));
}
if (i < 3)
{
IOUtil_1DLines (cctkGH, i + 1, origin_index, origin_phys, myGH->spxyz[i]);
}
}
/* free allocated resources */
free (origin_phys[0]);
free (origin_index[0]);
}
/*@@
@routine IOASCII_Choose2D
@author Gabrielle Allen
@date July 6 2000
@desc
Use parameters to choose the 2D slices through the output data.
@enddesc
@calls IOUtil_2DPlanes
@var GH
@vdesc Pointer to CCTK grid hierarchy
@vtype const cGH *
@vio in
@endvar
@@*/
void IOASCII_Choose2D (CCTK_ARGUMENTS)
{
int i, maxdim;
asciiioGH *myGH;
int origin_index[3];
CCTK_REAL origin_phys[3];
DECLARE_CCTK_PARAMETERS
GET_SLICE (out2D_xyplane_z, out_xyplane_z, origin_index[0], origin_phys[0]);
GET_SLICE (out2D_xzplane_y, out_xzplane_y, origin_index[1], origin_phys[1]);
GET_SLICE (out2D_yzplane_x, out_yzplane_x, origin_index[2], origin_phys[2]);
maxdim = CCTK_MaxDim ();
myGH = CCTK_GHExtension (cctkGH, "IOASCII");
myGH->sp2xyz = malloc (maxdim * sizeof (int *));
for (i = 0; i < maxdim; i++)
{
myGH->sp2xyz[i] = calloc (i + 1, sizeof (int));
if (i == 2)
{
IOUtil_2DPlanes (cctkGH, i + 1, origin_index, origin_phys, myGH->sp2xyz[i]);
}
}
}
|