| Geant4 Cross Reference |
1 /* 1 /*
2 # <<BEGIN-copyright>> 2 # <<BEGIN-copyright>>
3 # <<END-copyright>> 3 # <<END-copyright>>
4 */ 4 */
5 5
6 #include <map> 6 #include <map>
7 #include <string.h> 7 #include <string.h>
8 #include <cmath> 8 #include <cmath>
9 9
10 #include "MCGIDI.h" 10 #include "MCGIDI.h"
11 #include "MCGIDI_misc.h" 11 #include "MCGIDI_misc.h"
12 #include <xDataTOM_importXML_private.h> 12 #include <xDataTOM_importXML_private.h>
13 13
14 #if defined __cplusplus 14 #if defined __cplusplus
15 extern "C" { 15 extern "C" {
16 namespace GIDI { 16 namespace GIDI {
17 using namespace GIDI; 17 using namespace GIDI;
18 #endif 18 #endif
19 static int _MCGIDI_target_releaseAndReturnOne( 19 static int _MCGIDI_target_releaseAndReturnOne( statusMessageReporting *smr, MCGIDI_target *target );
20 #if defined __cplusplus 20 #if defined __cplusplus
21 } 21 }
22 } 22 }
23 #endif 23 #endif
24 /* 24 /*
25 ********************************************** 25 ************************************************************
26 */ 26 */
27 27
28 #if defined __cplusplus 28 #if defined __cplusplus
29 namespace GIDI { 29 namespace GIDI {
30 using namespace GIDI; 30 using namespace GIDI;
31 #endif 31 #endif
32 32
33 MCGIDI_target *MCGIDI_target_new( statusMessag 33 MCGIDI_target *MCGIDI_target_new( statusMessageReporting *smr ) {
34 34
35 MCGIDI_target *target; 35 MCGIDI_target *target;
36 36
37 if( ( target = (MCGIDI_target *) smr_mallo 37 if( ( target = (MCGIDI_target *) smr_malloc2( smr, sizeof( MCGIDI_target ), 0, "target" ) ) == NULL ) return( NULL );
38 if( MCGIDI_target_initialize( smr, target 38 if( MCGIDI_target_initialize( smr, target ) ) target = MCGIDI_target_free( smr, target );
39 return( target ); 39 return( target );
40 } 40 }
41 /* 41 /*
42 ********************************************** 42 ************************************************************
43 */ 43 */
44 int MCGIDI_target_initialize( statusMessageRep 44 int MCGIDI_target_initialize( statusMessageReporting * /*smr*/, MCGIDI_target *target ) {
45 45
46 memset( target, 0, sizeof( MCGIDI_target ) 46 memset( target, 0, sizeof( MCGIDI_target ) );
47 return( 0 ); 47 return( 0 );
48 } 48 }
49 /* 49 /*
50 ********************************************** 50 ************************************************************
51 */ 51 */
52 MCGIDI_target *MCGIDI_target_newRead( statusMe 52 MCGIDI_target *MCGIDI_target_newRead( statusMessageReporting *smr, const char *fileName ) {
53 53
54 MCGIDI_target *target; 54 MCGIDI_target *target;
55 55
56 if( ( target = MCGIDI_target_new( smr ) ) 56 if( ( target = MCGIDI_target_new( smr ) ) == NULL ) return( NULL );
57 if( MCGIDI_target_read( smr, target, fileN 57 if( MCGIDI_target_read( smr, target, fileName ) != 0 ) smr_freeMemory( (void **) &target );
58 return( target ); 58 return( target );
59 } 59 }
60 /* 60 /*
61 ********************************************** 61 ************************************************************
62 */ 62 */
63 int MCGIDI_target_readFromMapViaPoPIDs( status 63 int MCGIDI_target_readFromMapViaPoPIDs( statusMessageReporting *smr, MCGIDI_target *target, MCGIDI_map *map, const char *evaluation,
64 int projectile_PoPID, int target_PoPID 64 int projectile_PoPID, int target_PoPID ) {
65 65
66 char *targetPath; 66 char *targetPath;
67 67
68 if( ( targetPath = MCGIDI_map_findTargetVi 68 if( ( targetPath = MCGIDI_map_findTargetViaPoPIDs( smr, map, evaluation, projectile_PoPID, target_PoPID ) ) == NULL ) return( 1 );
69 return( MCGIDI_target_read( smr, target, t 69 return( MCGIDI_target_read( smr, target, targetPath ) );
70 } 70 }
71 /* 71 /*
72 ********************************************** 72 ************************************************************
73 */ 73 */
74 int MCGIDI_target_readFromMap( statusMessageRe 74 int MCGIDI_target_readFromMap( statusMessageReporting *smr, MCGIDI_target *target, MCGIDI_map *map, const char *evaluation, const char *projectileName,
75 const char *targetName ) { 75 const char *targetName ) {
76 76
77 char *targetPath; 77 char *targetPath;
78 78
79 if( ( targetPath = MCGIDI_map_findTarget( 79 if( ( targetPath = MCGIDI_map_findTarget( smr, map, evaluation, projectileName, targetName ) ) == NULL ) return( 1 );
80 return( MCGIDI_target_read( smr, target, t 80 return( MCGIDI_target_read( smr, target, targetPath ) );
81 } 81 }
82 /* 82 /*
83 ********************************************** 83 ************************************************************
84 */ 84 */
85 MCGIDI_target *MCGIDI_target_newReadFromMapVia 85 MCGIDI_target *MCGIDI_target_newReadFromMapViaPoPIDs( statusMessageReporting *smr, MCGIDI_map *map, const char *evaluation,
86 int projectile_PoPID, int target_PoPID 86 int projectile_PoPID, int target_PoPID ) {
87 87
88 char *targetPath; 88 char *targetPath;
89 MCGIDI_target *target; 89 MCGIDI_target *target;
90 90
91 if( ( targetPath = MCGIDI_map_findTargetVi 91 if( ( targetPath = MCGIDI_map_findTargetViaPoPIDs( smr, map, evaluation, projectile_PoPID, target_PoPID ) ) == NULL ) return( NULL );
92 target = MCGIDI_target_newRead( smr, targe 92 target = MCGIDI_target_newRead( smr, targetPath );
93 smr_freeMemory( (void **) &targetPath ); 93 smr_freeMemory( (void **) &targetPath );
94 return( target ); 94 return( target );
95 } 95 }
96 /* 96 /*
97 ********************************************** 97 ************************************************************
98 */ 98 */
99 MCGIDI_target *MCGIDI_target_newReadFromMap( s 99 MCGIDI_target *MCGIDI_target_newReadFromMap( statusMessageReporting *smr, MCGIDI_map *map, const char *evaluation, const char *projectileName,
100 const char *targetName ) { 100 const char *targetName ) {
101 101
102 char *targetPath; 102 char *targetPath;
103 MCGIDI_target *target; 103 MCGIDI_target *target;
104 104
105 targetPath = MCGIDI_map_findTarget( smr, m 105 targetPath = MCGIDI_map_findTarget( smr, map, evaluation, projectileName, targetName );
106 if( targetPath == NULL ) return( NULL ); 106 if( targetPath == NULL ) return( NULL );
107 target = MCGIDI_target_newRead( smr, targe 107 target = MCGIDI_target_newRead( smr, targetPath );
108 smr_freeMemory( (void **) &targetPath ); 108 smr_freeMemory( (void **) &targetPath );
109 return( target ); 109 return( target );
110 } 110 }
111 /* 111 /*
112 ********************************************** 112 ************************************************************
113 */ 113 */
114 MCGIDI_target *MCGIDI_target_free( statusMessa 114 MCGIDI_target *MCGIDI_target_free( statusMessageReporting *smr, MCGIDI_target *target ) {
115 115
116 MCGIDI_target_release( smr, target ); 116 MCGIDI_target_release( smr, target );
117 smr_freeMemory( (void **) &target ); 117 smr_freeMemory( (void **) &target );
118 return( NULL ); 118 return( NULL );
119 } 119 }
120 /* 120 /*
121 ********************************************** 121 ************************************************************
122 */ 122 */
123 int MCGIDI_target_release( statusMessageReport 123 int MCGIDI_target_release( statusMessageReporting *smr, MCGIDI_target *target ) {
124 124
125 int i; 125 int i;
126 126
127 smr_freeMemory( (void **) &(target->path) 127 smr_freeMemory( (void **) &(target->path) );
128 smr_freeMemory( (void **) &(target->absPat 128 smr_freeMemory( (void **) &(target->absPath) );
129 xDataTOMAL_release( &(target->attributes) 129 xDataTOMAL_release( &(target->attributes) );
130 for( i = 0; i < target->nHeatedTargets; i+ 130 for( i = 0; i < target->nHeatedTargets; i++ ) {
131 smr_freeMemory( (void **) &(target->he 131 smr_freeMemory( (void **) &(target->heatedTargets[i].path) );
132 smr_freeMemory( (void **) &(target->he 132 smr_freeMemory( (void **) &(target->heatedTargets[i].contents) );
133 if( target->heatedTargets[i].heatedTar 133 if( target->heatedTargets[i].heatedTarget != NULL ) MCGIDI_target_heated_free( smr, target->heatedTargets[i].heatedTarget );
134 } 134 }
135 smr_freeMemory( (void **) &(target->heated 135 smr_freeMemory( (void **) &(target->heatedTargets) );
136 smr_freeMemory( (void **) &(target->readHe 136 smr_freeMemory( (void **) &(target->readHeatedTargets) );
137 MCGIDI_target_initialize( smr, target ); 137 MCGIDI_target_initialize( smr, target );
138 return( 0 ); 138 return( 0 );
139 } 139 }
140 /* 140 /*
141 ********************************************** 141 ************************************************************
142 */ 142 */
143 int MCGIDI_target_read( statusMessageReporting 143 int MCGIDI_target_read( statusMessageReporting *smr, MCGIDI_target *target, const char *fileName ) {
144 /* 144 /*
145 * If a target has already been read into thi 145 * If a target has already been read into this target, user must have called MCGIDI_target_release before calling this routine.
146 * Otherwise, there will be memory leaks. 146 * Otherwise, there will be memory leaks.
147 */ 147 */
148 xDataXML_document *doc; 148 xDataXML_document *doc;
149 xDataXML_element *element, *child; 149 xDataXML_element *element, *child;
150 int i, iHeated, nHeated = 0, status = 1; 150 int i, iHeated, nHeated = 0, status = 1;
151 double temperature; 151 double temperature;
152 /* char *pReturnValue; */ 152 /* char *pReturnValue; */
153 char const *version, *contents; 153 char const *version, *contents;
154 154
155 MCGIDI_target_initialize( smr, target ); 155 MCGIDI_target_initialize( smr, target );
156 if( ( target->path = smr_allocateCopyStrin 156 if( ( target->path = smr_allocateCopyString2( smr, fileName, "path" ) ) == NULL ) return( status );
157 if( ( target->absPath = MCGIDI_misc_getAbs 157 if( ( target->absPath = MCGIDI_misc_getAbsPath( smr, fileName ) ) == NULL ) return( _MCGIDI_target_releaseAndReturnOne( smr, target ) );
158 if( ( doc = xDataXML_importFile2( smr, fil 158 if( ( doc = xDataXML_importFile2( smr, fileName ) ) == NULL ) return( _MCGIDI_target_releaseAndReturnOne( smr, target ) );
159 element = xDataXML_getDocumentsElement( do 159 element = xDataXML_getDocumentsElement( doc );
160 if( strcmp( element->name, "xTarget" ) != 160 if( strcmp( element->name, "xTarget" ) != 0 ) {
161 MCGIDI_misc_setMessageError_Element( s 161 MCGIDI_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "input file's top element must be xTarget and not %s", element->name ); }
162 else { 162 else {
163 status = 0; 163 status = 0;
164 if( ( version = xDataXML_getAttributes 164 if( ( version = xDataXML_getAttributesValueInElement( element, "version" ) ) == NULL ) {
165 smr_setReportError2( smr, smr_unkn 165 smr_setReportError2( smr, smr_unknownID, 1, "version attribute missing from element '%s'", element->name );
166 status = 1; } 166 status = 1; }
167 else { 167 else {
168 if( strcmp( version, "xMCProcess 0 168 if( strcmp( version, "xMCProcess 0.1" ) != 0 ) {
169 smr_setReportError2( smr, smr_ 169 smr_setReportError2( smr, smr_unknownID, 1, "Unsupported version '%s' for element %s", version, element->name );
170 status = 1; 170 status = 1;
171 } 171 }
172 } 172 }
173 if( status == 0 ) { 173 if( status == 0 ) {
174 /* pReturnValue = ( MCGIDI_misc_cop 174 /* pReturnValue = ( MCGIDI_misc_copyXMLAttributesToTOM( smr, &(target->attributes), &(element->attributes) ) ) ? NULL : target->path; */
175 for( nHeated = 0, child = xDataXML 175 for( nHeated = 0, child = xDataXML_getFirstElement( element ); child != NULL; nHeated++, child = xDataXML_getNextElement( child ) ) {
176 if( strcmp( child->name, "targ 176 if( strcmp( child->name, "target" ) != 0 ) {
177 MCGIDI_misc_setMessageErro 177 MCGIDI_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "element can only have target sub-elements%s",
178 element->name 178 element->name );
179 status = 1; 179 status = 1;
180 break; 180 break;
181 } 181 }
182 } 182 }
183 } 183 }
184 if( status == 0 ) { 184 if( status == 0 ) {
185 if( ( target->heatedTargets = (MCG 185 if( ( target->heatedTargets = (MCGIDI_target_heated_info *) smr_malloc2( smr, nHeated * sizeof( MCGIDI_target_heated_info ), 1, "heatedTargets" ) ) == NULL ) {
186 status = 1; } 186 status = 1; }
187 else { 187 else {
188 if( ( target->readHeatedTarget 188 if( ( target->readHeatedTargets = (MCGIDI_target_heated_info **) smr_malloc2( smr, nHeated * sizeof( MCGIDI_target_heated_info * ), 1, "heatedTargets" ) ) == NULL )
189 status = 1; 189 status = 1;
190 } 190 }
191 for( nHeated = 0, child = xDataXML 191 for( nHeated = 0, child = xDataXML_getFirstElement( element ); ( status == 0 ) && ( child != NULL );
192 nHeated++, child = xDataXM 192 nHeated++, child = xDataXML_getNextElement( child ) ) {
193 if( ( i = xDataXML_convertAttr 193 if( ( i = xDataXML_convertAttributeToDouble( smr, child, "temperature", &temperature, 1 ) ) != 0 ) {
194 if( i > 0 ) smr_setReportE 194 if( i > 0 ) smr_setReportError2p( smr, smr_unknownID, 1, "target does not have a temperature attribute" );
195 status = 1; 195 status = 1;
196 break; 196 break;
197 } 197 }
198 for( iHeated = 0; iHeated < nH 198 for( iHeated = 0; iHeated < nHeated; iHeated++ ) if( target->heatedTargets[iHeated].temperature > temperature ) break;
199 if( iHeated < nHeated ) for( i 199 if( iHeated < nHeated ) for( i = nHeated; i >= iHeated; i-- ) target->heatedTargets[i+1] = target->heatedTargets[i];
200 target->heatedTargets[iHeated] 200 target->heatedTargets[iHeated].temperature = temperature;
201 target->heatedTargets[iHeated] 201 target->heatedTargets[iHeated].path = NULL;
202 target->heatedTargets[iHeated] 202 target->heatedTargets[iHeated].contents = NULL;
203 target->heatedTargets[iHeated] 203 target->heatedTargets[iHeated].heatedTarget = NULL;
204 if( ( contents = xDataXML_getA 204 if( ( contents = xDataXML_getAttributesValueInElement( child, "contents" ) ) != NULL ) {
205 if( ( target->heatedTarget 205 if( ( target->heatedTargets[iHeated].contents = smr_allocateCopyString2( smr, contents, "contents" ) ) == NULL ) {
206 status = 1; 206 status = 1;
207 break; 207 break;
208 } 208 }
209 } 209 }
210 if( ( contents = xDataXML_getA 210 if( ( contents = xDataXML_getAttributesValueInElement( child, "file" ) ) == NULL ) {
211 status = 1; 211 status = 1;
212 break; 212 break;
213 } 213 }
214 if( ( target->heatedTargets[iH 214 if( ( target->heatedTargets[iHeated].path = (char *) smr_malloc2(smr, strlen( target->absPath ) + strlen( contents ) + 2, 0, "path") ) == NULL) {
215 status = 1; 215 status = 1;
216 break; 216 break;
217 } 217 }
218 strcpy( target->heatedTargets[ 218 strcpy( target->heatedTargets[iHeated].path, target->absPath );
219 *strrchr( target->heatedTarget 219 *strrchr( target->heatedTargets[iHeated].path, '/' ) = 0;
220 strcat( target->heatedTargets[ 220 strcat( target->heatedTargets[iHeated].path, "/" );
221 strcat( target->heatedTargets[ 221 strcat( target->heatedTargets[iHeated].path, contents );
222 target->nHeatedTargets++; 222 target->nHeatedTargets++;
223 } 223 }
224 } 224 }
225 } 225 }
226 xDataXML_freeDoc( smr, doc ); 226 xDataXML_freeDoc( smr, doc );
227 if( status == 0 ) { 227 if( status == 0 ) {
228 for( i = 0; i < nHeated; i++ ) target- 228 for( i = 0; i < nHeated; i++ ) target->heatedTargets[i].ordinal = i;
229 for( i = 0; i < nHeated; i++ ) if( tar 229 for( i = 0; i < nHeated; i++ ) if( target->heatedTargets[i].contents == NULL ) break;
230 if( i == nHeated ) i = 0; 230 if( i == nHeated ) i = 0; /* All heated targets are crossSection only. */
231 if( MCGIDI_target_readHeatedTarget( sm 231 if( MCGIDI_target_readHeatedTarget( smr, target, i ) == 0 ) {
232 target->baseHeatedTarget = target- 232 target->baseHeatedTarget = target->heatedTargets[i].heatedTarget; }
233 else { 233 else {
234 MCGIDI_target_release( NULL, targe 234 MCGIDI_target_release( NULL, target );
235 status = 1; 235 status = 1;
236 } } 236 } }
237 else { 237 else {
238 MCGIDI_target_release( smr, target ); 238 MCGIDI_target_release( smr, target );
239 } 239 }
240 return( status ); 240 return( status );
241 } 241 }
242 /* 242 /*
243 ********************************************** 243 ************************************************************
244 */ 244 */
245 char const *MCGIDI_target_getAttributesValue( 245 char const *MCGIDI_target_getAttributesValue( statusMessageReporting * /*smr*/, MCGIDI_target *target, char const *name ) {
246 246
247 return( xDataTOMAL_getAttributesValue( &(t 247 return( xDataTOMAL_getAttributesValue( &(target->attributes), name ) );
248 } 248 }
249 /* 249 /*
250 ********************************************** 250 ************************************************************
251 */ 251 */
252 int MCGIDI_target_getTemperatures( statusMessa 252 int MCGIDI_target_getTemperatures( statusMessageReporting * /*smr*/, MCGIDI_target *target, double *temperatures ) {
253 253
254 int i; 254 int i;
255 255
256 if( temperatures != NULL ) for( i = 0; i < 256 if( temperatures != NULL ) for( i = 0; i < target->nHeatedTargets; i++ ) temperatures[i] = target->heatedTargets[i].temperature;
257 return( target->nHeatedTargets ); 257 return( target->nHeatedTargets );
258 } 258 }
259 /* 259 /*
260 ********************************************** 260 ************************************************************
261 */ 261 */
262 int MCGIDI_target_readHeatedTarget( statusMess 262 int MCGIDI_target_readHeatedTarget( statusMessageReporting *smr, MCGIDI_target *target, int index ) {
263 263
264 int i; 264 int i;
265 265
266 if( ( index < 0 ) || ( index >= target->nH 266 if( ( index < 0 ) || ( index >= target->nHeatedTargets ) ) {
267 smr_setReportError2( smr, smr_unknownI 267 smr_setReportError2( smr, smr_unknownID, 1, "temperature index = %d out of range (0 <= index < %d", index, target->nHeatedTargets );
268 return( -1 ); 268 return( -1 );
269 } 269 }
270 if( target->heatedTargets[index].heatedTar 270 if( target->heatedTargets[index].heatedTarget != NULL ) return( 1 );
271 if( ( target->heatedTargets[index].heatedT 271 if( ( target->heatedTargets[index].heatedTarget = MCGIDI_target_heated_newRead( smr, target->heatedTargets[index].path ) ) != NULL ) {
272 target->projectilePOP = target->heated 272 target->projectilePOP = target->heatedTargets[index].heatedTarget->projectilePOP;
273 target->targetPOP = target->heatedTarg 273 target->targetPOP = target->heatedTargets[index].heatedTarget->targetPOP;
274 if( target->heatedTargets[index].heate 274 if( target->heatedTargets[index].heatedTarget != NULL ) {
275 target->heatedTargets[index].heate 275 target->heatedTargets[index].heatedTarget->ordinal = target->heatedTargets[index].ordinal;
276 for( i = target->nReadHeatedTarget 276 for( i = target->nReadHeatedTargets; i > 0; i-- ) {
277 if( target->readHeatedTargets[ 277 if( target->readHeatedTargets[i-1]->temperature < target->heatedTargets[index].temperature ) break;
278 target->readHeatedTargets[i] = 278 target->readHeatedTargets[i] = target->readHeatedTargets[i-1];
279 } 279 }
280 target->readHeatedTargets[i] = &(t 280 target->readHeatedTargets[i] = &(target->heatedTargets[i]);
281 target->nReadHeatedTargets++; 281 target->nReadHeatedTargets++;
282 } 282 }
283 } 283 }
284 return( ( target->heatedTargets[index].hea 284 return( ( target->heatedTargets[index].heatedTarget == NULL ? -1 : 0 ) );
285 } 285 }
286 /* 286 /*
287 ********************************************** 287 ************************************************************
288 */ 288 */
289 MCGIDI_target_heated *MCGIDI_target_getHeatedT 289 MCGIDI_target_heated *MCGIDI_target_getHeatedTargetAtIndex_ReadIfNeeded( statusMessageReporting *smr, MCGIDI_target *target, int index ) {
290 290
291 if( ( index < 0 ) || ( index >= target->nH 291 if( ( index < 0 ) || ( index >= target->nHeatedTargets ) ) {
292 smr_setReportError2( smr, smr_unknownI 292 smr_setReportError2( smr, smr_unknownID, 1, "temperature index = %d out of range (0 <= index < %d", index, target->nHeatedTargets );
293 return( NULL ); 293 return( NULL );
294 } 294 }
295 if( target->heatedTargets[index].heatedTar 295 if( target->heatedTargets[index].heatedTarget == NULL ) MCGIDI_target_readHeatedTarget( smr, target, index );
296 return( target->heatedTargets[index].heate 296 return( target->heatedTargets[index].heatedTarget );
297 } 297 }
298 /* 298 /*
299 ********************************************** 299 ************************************************************
300 */ 300 */
301 MCGIDI_target_heated *MCGIDI_target_getHeatedT 301 MCGIDI_target_heated *MCGIDI_target_getHeatedTargetAtTIndex( statusMessageReporting *smr, MCGIDI_target *target, int index ) {
302 302
303 if( ( index < 0 ) || ( index >= target->nH 303 if( ( index < 0 ) || ( index >= target->nHeatedTargets ) ) {
304 smr_setReportError2( smr, smr_unknownI 304 smr_setReportError2( smr, smr_unknownID, 1, "temperature index = %d out of range (0 <= index < %d", index, target->nHeatedTargets );
305 return( NULL ); 305 return( NULL );
306 } 306 }
307 if( target->heatedTargets[index].heatedTar 307 if( target->heatedTargets[index].heatedTarget == NULL ) {
308 smr_setReportError2( smr, smr_unknownI 308 smr_setReportError2( smr, smr_unknownID, 1, "temperature index = %d not read in", index );
309 return( NULL ); 309 return( NULL );
310 } 310 }
311 return( target->heatedTargets[index].heate 311 return( target->heatedTargets[index].heatedTarget );
312 } 312 }
313 /* 313 /*
314 ********************************************** 314 ************************************************************
315 */ 315 */
316 int MCGIDI_target_numberOfReactions( statusMes 316 int MCGIDI_target_numberOfReactions( statusMessageReporting *smr, MCGIDI_target *target ) {
317 317
318 return( MCGIDI_target_heated_numberOfReact 318 return( MCGIDI_target_heated_numberOfReactions( smr, target->baseHeatedTarget ) );
319 } 319 }
320 /* 320 /*
321 ********************************************** 321 ************************************************************
322 */ 322 */
323 enum MCGIDI_reactionType MCGIDI_target_getReac 323 enum MCGIDI_reactionType MCGIDI_target_getReactionTypeAtIndex( statusMessageReporting *smr, MCGIDI_target *target, int index ) {
324 324
325 MCGIDI_reaction *reaction = MCGIDI_target_ 325 MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex_smr( smr, target->baseHeatedTarget, index );
326 326
327 if( reaction == NULL ) return( MCGIDI_reac 327 if( reaction == NULL ) return( MCGIDI_reactionType_unknown_e );
328 return( MCGIDI_reaction_getReactionType( s 328 return( MCGIDI_reaction_getReactionType( smr, reaction ) );
329 } 329 }
330 /* 330 /*
331 ********************************************** 331 ************************************************************
332 */ 332 */
333 MCGIDI_reaction *MCGIDI_target_getReactionAtIn 333 MCGIDI_reaction *MCGIDI_target_getReactionAtIndex( MCGIDI_target *target, int index ) {
334 334
335 return( MCGIDI_target_heated_getReactionAt 335 return( MCGIDI_target_heated_getReactionAtIndex( target->baseHeatedTarget, index ) );
336 } 336 }
337 /* 337 /*
338 ********************************************** 338 ************************************************************
339 */ 339 */
340 MCGIDI_reaction *MCGIDI_target_getReactionAtIn 340 MCGIDI_reaction *MCGIDI_target_getReactionAtIndex_smr( statusMessageReporting *smr, MCGIDI_target *target, int index ) {
341 341
342 return( MCGIDI_target_heated_getReactionAt 342 return( MCGIDI_target_heated_getReactionAtIndex_smr( smr, target->baseHeatedTarget, index ) );
343 } 343 }
344 /* 344 /*
345 ********************************************** 345 ************************************************************
346 */ 346 */
347 int MCGIDI_target_numberOfProductionReactions( 347 int MCGIDI_target_numberOfProductionReactions( statusMessageReporting * /*smr*/, MCGIDI_target * /*target*/ ) {
348 348
349 #if 0 349 #if 0
350 return( MCGIDI_target_heated_numberOfProdu 350 return( MCGIDI_target_heated_numberOfProductionReactions( smr, target->baseHeatedTarget ) );
351 #endif 351 #endif
352 return( 0 ); 352 return( 0 );
353 } 353 }
354 354
355 /* 355 /*
356 ********************************************** 356 ************************************************************
357 */ 357 */
358 double MCGIDI_target_getTotalCrossSectionAtTAn 358 double MCGIDI_target_getTotalCrossSectionAtTAndE( statusMessageReporting *smr, MCGIDI_target *target, MCGIDI_quantitiesLookupModes &modes,
359 bool sampling ) { 359 bool sampling ) {
360 360
361 int i; 361 int i;
362 double xsec = 0., xsec1, xsec2, temperatur 362 double xsec = 0., xsec1, xsec2, temperature = modes.getTemperature( );
363 363
364 for( i = 0; i < target->nReadHeatedTargets 364 for( i = 0; i < target->nReadHeatedTargets; i++ ) if( target->readHeatedTargets[i]->temperature > temperature ) break;
365 if( i == 0 ) { 365 if( i == 0 ) {
366 xsec = MCGIDI_target_heated_getTotalCr 366 xsec = MCGIDI_target_heated_getTotalCrossSectionAtE( smr, target->readHeatedTargets[0]->heatedTarget, modes, sampling ); }
367 else if( i == target->nReadHeatedTargets ) 367 else if( i == target->nReadHeatedTargets ) {
368 xsec = MCGIDI_target_heated_getTotalCr 368 xsec = MCGIDI_target_heated_getTotalCrossSectionAtE( smr, target->readHeatedTargets[i-1]->heatedTarget, modes, sampling ); }
369 else { 369 else {
370 xsec1 = MCGIDI_target_heated_getTotalC 370 xsec1 = MCGIDI_target_heated_getTotalCrossSectionAtE( smr, target->readHeatedTargets[i-1]->heatedTarget, modes, sampling );
371 xsec2 = MCGIDI_target_heated_getTotalC 371 xsec2 = MCGIDI_target_heated_getTotalCrossSectionAtE( smr, target->readHeatedTargets[i ]->heatedTarget, modes, sampling );
372 xsec = ( ( target->readHeatedTargets[i 372 xsec = ( ( target->readHeatedTargets[i]->temperature - temperature ) * xsec1 +
373 ( temperature - target->readH 373 ( temperature - target->readHeatedTargets[i-1]->temperature ) * xsec2 ) /
374 ( target->readHeatedTargets[i 374 ( target->readHeatedTargets[i]->temperature - target->readHeatedTargets[i-1]->temperature );
375 } 375 }
376 376
377 return( xsec ); 377 return( xsec );
378 } 378 }
379 /* 379 /*
380 ********************************************** 380 ************************************************************
381 */ 381 */
382 int MCGIDI_target_getDomain( statusMessageRepo 382 int MCGIDI_target_getDomain( statusMessageReporting *smr, MCGIDI_target *target, double *EMin, double *EMax ) {
383 383
384 int ir, nr = MCGIDI_target_numberOfReactio 384 int ir, nr = MCGIDI_target_numberOfReactions( smr, target );
385 double EMin_, EMax_; 385 double EMin_, EMax_;
386 386
387 for( ir = 0; ir < nr; ir++ ) { 387 for( ir = 0; ir < nr; ir++ ) {
388 MCGIDI_target_heated_getReactionsDomai 388 MCGIDI_target_heated_getReactionsDomain( smr, target->baseHeatedTarget, ir, &EMin_, &EMax_ );
389 if( ir == 0 ) { 389 if( ir == 0 ) {
390 *EMin = EMin_; 390 *EMin = EMin_;
391 *EMax = EMax_; } 391 *EMax = EMax_; }
392 else { 392 else {
393 if( *EMin > EMin_ ) *EMin = EMin_; 393 if( *EMin > EMin_ ) *EMin = EMin_;
394 if( *EMax < EMax_ ) *EMax = EMax_; 394 if( *EMax < EMax_ ) *EMax = EMax_;
395 } 395 }
396 } 396 }
397 return( 0 ); 397 return( 0 );
398 } 398 }
399 /* 399 /*
400 ********************************************** 400 ************************************************************
401 */ 401 */
402 double MCGIDI_target_getIndexReactionCrossSect 402 double MCGIDI_target_getIndexReactionCrossSectionAtE( statusMessageReporting *smr, MCGIDI_target *target, int index, MCGIDI_quantitiesLookupModes &modes,
403 bool sampling ) { 403 bool sampling ) {
404 404
405 int i; 405 int i;
406 double xsec = 0., xsec1, xsec2, temperatur 406 double xsec = 0., xsec1, xsec2, temperature = modes.getTemperature( );
407 407
408 for( i = 0; i < target->nReadHeatedTargets 408 for( i = 0; i < target->nReadHeatedTargets; i++ ) if( target->readHeatedTargets[i]->temperature > temperature ) break;
409 if( i == 0 ) { 409 if( i == 0 ) {
410 xsec = MCGIDI_target_heated_getIndexRe 410 xsec = MCGIDI_target_heated_getIndexReactionCrossSectionAtE( smr, target->readHeatedTargets[0]->heatedTarget, index, modes, sampling ); }
411 else if( i == target->nReadHeatedTargets ) 411 else if( i == target->nReadHeatedTargets ) {
412 xsec = MCGIDI_target_heated_getIndexRe 412 xsec = MCGIDI_target_heated_getIndexReactionCrossSectionAtE( smr, target->readHeatedTargets[i-1]->heatedTarget, index, modes, sampling ); }
413 else { 413 else {
414 xsec1 = MCGIDI_target_heated_getIndexR 414 xsec1 = MCGIDI_target_heated_getIndexReactionCrossSectionAtE(smr, target->readHeatedTargets[i-1]->heatedTarget, index, modes, sampling );
415 xsec2 = MCGIDI_target_heated_getIndexR 415 xsec2 = MCGIDI_target_heated_getIndexReactionCrossSectionAtE(smr, target->readHeatedTargets[i ]->heatedTarget, index, modes, sampling );
416 xsec = ( ( target->readHeatedTargets[i 416 xsec = ( ( target->readHeatedTargets[i]->temperature - temperature ) * xsec1 +
417 ( temperature - target->readH 417 ( temperature - target->readHeatedTargets[i-1]->temperature ) * xsec2 ) /
418 ( target->readHeatedTargets[i 418 ( target->readHeatedTargets[i]->temperature - target->readHeatedTargets[i-1]->temperature );
419 } 419 }
420 420
421 return( xsec ); 421 return( xsec );
422 } 422 }
423 /* 423 /*
424 ********************************************** 424 ************************************************************
425 */ 425 */
426 int MCGIDI_target_sampleReaction( statusMessag 426 int MCGIDI_target_sampleReaction( statusMessageReporting *smr, MCGIDI_target *target, MCGIDI_quantitiesLookupModes &modes, double totalXSec,
427 double (*userrng)( void * ), void *rng 427 double (*userrng)( void * ), void *rngState ) {
428 428
429 int ir, nr = MCGIDI_target_numberOfReactio 429 int ir, nr = MCGIDI_target_numberOfReactions( smr, target );
430 double rngValue = (*userrng)( rngState ); 430 double rngValue = (*userrng)( rngState );
431 double cumm_xsec = 0., r_xsec = rngValue * 431 double cumm_xsec = 0., r_xsec = rngValue * totalXSec;
432 432
433 for( ir = 0; ir < nr; ir++ ) { 433 for( ir = 0; ir < nr; ir++ ) {
434 cumm_xsec += MCGIDI_target_getIndexRea 434 cumm_xsec += MCGIDI_target_getIndexReactionCrossSectionAtE( smr, target, ir, modes, true );
435 if( cumm_xsec >= r_xsec ) break; 435 if( cumm_xsec >= r_xsec ) break;
436 } 436 }
437 if( ir == nr ) { 437 if( ir == nr ) {
438 if( ( totalXSec - cumm_xsec ) >= 1e-12 438 if( ( totalXSec - cumm_xsec ) >= 1e-12 * totalXSec ) {
439 smr_setReportError2( smr, smr_unkn 439 smr_setReportError2( smr, smr_unknownID, 1,
440 "Failed to sample a reaction f 440 "Failed to sample a reaction for temperature = %.12e, energy = %.12e, totalXSec = %16.e, rngValue = %16.e, r_xsec = %16.e, cumm_xsec = %16.e",
441 modes.getTemperature( ), modes 441 modes.getTemperature( ), modes.getProjectileEnergy( ), totalXSec, rngValue, r_xsec, cumm_xsec );
442 return( -1 ); 442 return( -1 );
443 } 443 }
444 ir--; /* 444 ir--; /* May not be correct but close. */
445 } 445 }
446 if( modes.getCrossSectionMode( ) == MCGIDI 446 if( modes.getCrossSectionMode( ) == MCGIDI_quantityLookupMode_grouped ) {
447 MCGIDI_reaction *reaction = MCGIDI_tar 447 MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex( target->baseHeatedTarget, ir );
448 448
449 if( modes.getGroupIndex( ) == reaction 449 if( modes.getGroupIndex( ) == reaction->thresholdGroupIndex ) {
450 double dEnergy = modes.getProjecti 450 double dEnergy = modes.getProjectileEnergy( ) - reaction->EMin;
451 451
452 if( dEnergy <= 0 ) return( MCGIDI 452 if( dEnergy <= 0 ) return( MCGIDI_nullReaction );
453 if( ( (*userrng)( rngState ) * rea 453 if( ( (*userrng)( rngState ) * reaction->thresholdGroupDomain ) > dEnergy ) return( MCGIDI_nullReaction );
454 } 454 }
455 } 455 }
456 return( ir ); 456 return( ir );
457 } 457 }
458 /* 458 /*
459 ********************************************** 459 ************************************************************
460 */ 460 */
461 int MCGIDI_target_sampleNullReactionProductsAt 461 int MCGIDI_target_sampleNullReactionProductsAtE( statusMessageReporting *smr, MCGIDI_target *target,
462 MCGIDI_quantitiesLookupModes &modes, MCGID 462 MCGIDI_quantitiesLookupModes &modes, MCGIDI_decaySamplingInfo *decaySamplingInfo, MCGIDI_sampledProductsDatas *productDatas ) {
463 463
464 MCGIDI_sampledProductsData productData; 464 MCGIDI_sampledProductsData productData;
465 465
466 productData.isVelocity = decaySamplingInfo 466 productData.isVelocity = decaySamplingInfo->isVelocity;
467 productData.pop = target->projectilePOP; 467 productData.pop = target->projectilePOP;
468 productData.kineticEnergy = modes.getProje 468 productData.kineticEnergy = modes.getProjectileEnergy( );
469 productData.px_vx = 0.; 469 productData.px_vx = 0.;
470 productData.py_vy = 0.; 470 productData.py_vy = 0.;
471 productData.pz_vz = std::sqrt( productData 471 productData.pz_vz = std::sqrt( productData.kineticEnergy * ( productData.kineticEnergy + 2. * productData.pop->mass_MeV ) );
472 if( productData.isVelocity ) productData.p 472 if( productData.isVelocity ) productData.pz_vz *=
473 MCGIDI_speedOfLight_cm_sec / std:: 473 MCGIDI_speedOfLight_cm_sec / std::sqrt( productData.pz_vz * productData.pz_vz + productData.pop->mass_MeV * productData.pop->mass_MeV );
474 productData.delayedNeutronIndex = 0; 474 productData.delayedNeutronIndex = 0;
475 productData.delayedNeutronRate = 0.; 475 productData.delayedNeutronRate = 0.;
476 productData.birthTimeSec = 0; 476 productData.birthTimeSec = 0;
477 477
478 productDatas->numberOfProducts = 0; 478 productDatas->numberOfProducts = 0;
479 MCGIDI_sampledProducts_addProduct( smr, pr 479 MCGIDI_sampledProducts_addProduct( smr, productDatas, &productData );
480 return( productDatas->numberOfProducts ); 480 return( productDatas->numberOfProducts );
481 } 481 }
482 /* 482 /*
483 ********************************************** 483 ************************************************************
484 */ 484 */
485 int MCGIDI_target_sampleIndexReactionProductsA 485 int MCGIDI_target_sampleIndexReactionProductsAtE( statusMessageReporting *smr, MCGIDI_target *target, int index,
486 MCGIDI_quantitiesLookupModes &modes, M 486 MCGIDI_quantitiesLookupModes &modes, MCGIDI_decaySamplingInfo *decaySamplingInfo, MCGIDI_sampledProductsDatas *productData ) {
487 487
488 return( MCGIDI_target_heated_sampleIndexRe 488 return( MCGIDI_target_heated_sampleIndexReactionProductsAtE( smr, target->baseHeatedTarget, index, modes, decaySamplingInfo, productData ) );
489 } 489 }
490 /* 490 /*
491 ********************************************** 491 ************************************************************
492 */ 492 */
493 double MCGIDI_target_getIndexReactionFinalQ( s 493 double MCGIDI_target_getIndexReactionFinalQ( statusMessageReporting *smr, MCGIDI_target *target, int index,
494 MCGIDI_quantitiesLookupModes &modes ) 494 MCGIDI_quantitiesLookupModes &modes ) {
495 495
496 return( MCGIDI_target_heated_getIndexReact 496 return( MCGIDI_target_heated_getIndexReactionFinalQ( smr, target->baseHeatedTarget, index, modes ) );
497 } 497 }
498 /* 498 /*
499 ********************************************** 499 ************************************************************
500 */ 500 */
501 std::map<int, enum MCGIDI_transportability> co 501 std::map<int, enum MCGIDI_transportability> const *MCGIDI_target_getUniqueProducts( statusMessageReporting *smr, MCGIDI_target *target ) {
502 502
503 return( MCGIDI_target_heated_getUniqueProd 503 return( MCGIDI_target_heated_getUniqueProducts( smr, target->baseHeatedTarget ) );
504 } 504 }
505 /* 505 /*
506 ********************************************** 506 ************************************************************
507 */ 507 */
508 int MCGIDI_target_recast( statusMessageReporti 508 int MCGIDI_target_recast( statusMessageReporting *smr, MCGIDI_target *target, GIDI_settings &settings ) {
509 509
510 int i1, status = 0; 510 int i1, status = 0;
511 511
512 for( i1 = 0; i1 < target->nReadHeatedTarge 512 for( i1 = 0; i1 < target->nReadHeatedTargets; i1++ ) {
513 if( ( status = MCGIDI_target_heated_re 513 if( ( status = MCGIDI_target_heated_recast( smr, target->readHeatedTargets[i1]->heatedTarget, settings ) ) != 0 ) break;
514 } 514 }
515 return( status ); 515 return( status );
516 } 516 }
517 /* 517 /*
518 ********************************************** 518 ************************************************************
519 */ 519 */
520 static int _MCGIDI_target_releaseAndReturnOne( 520 static int _MCGIDI_target_releaseAndReturnOne( statusMessageReporting *smr, MCGIDI_target *target ) {
521 521
522 MCGIDI_target_release( smr, target ); 522 MCGIDI_target_release( smr, target );
523 return( 1 ); 523 return( 1 );
524 } 524 }
525 525
526 #if defined __cplusplus 526 #if defined __cplusplus
527 } 527 }
528 #endif 528 #endif
529 529