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