Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/lend/src/MCGIDI_target_heated.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  1 /*
  2 # <<BEGIN-copyright>>
  3 # <<END-copyright>>
  4 */
  5 
  6 #include <map>
  7 #include <string.h>
  8 #include <cmath>
  9 
 10 #include <xDataTOM.h>
 11 #include "MCGIDI.h"
 12 #include "MCGIDI_misc.h"
 13 #include "MCGIDI_private.h"
 14 
 15 #if defined __cplusplus
 16 namespace GIDI {
 17 using namespace GIDI;
 18 #endif
 19 
 20 static int MCGIDI_target_heated_parsePOPs( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element,
 21     xDataTOM_element *particleAliases );
 22 static int MCGIDI_target_heated_parseParticle( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element,
 23     xDataTOM_element *particleAliases );
 24 static int MCGIDI_target_heated_parseParticleLevel( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element, MCGIDI_POP *parent,
 25     double mass_MeV, xDataTOM_element *particleAliases );
 26 static int MCGIDI_target_heated_parseParticleGammas( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element, char const *name );
 27 static int MCGIDI_target_heated_parseReaction( statusMessageReporting *smr, xDataTOM_element *child, MCGIDI_target_heated *target, 
 28     MCGIDI_POPs *pops, MCGIDI_reaction *reaction );
 29 /*
 30 ************************************************************
 31 */
 32 MCGIDI_target_heated *MCGIDI_target_heated_new( statusMessageReporting *smr ) {
 33 
 34     MCGIDI_target_heated *target;
 35 
 36     if( ( target = (MCGIDI_target_heated *) smr_malloc2( smr, sizeof( MCGIDI_target_heated ), 0, "target" ) ) == NULL ) return( NULL );
 37     if( MCGIDI_target_heated_initialize( smr, target ) ) target = (MCGIDI_target_heated *) smr_freeMemory( (void **) &target );
 38     return( target );
 39 }
 40 /*
 41 ************************************************************
 42 */
 43 int MCGIDI_target_heated_initialize( statusMessageReporting *smr, MCGIDI_target_heated *target ) {
 44 
 45     memset( target, 0, sizeof( MCGIDI_target_heated ) );
 46     MCGIDI_POPs_initial( smr, &(target->pops), 100 );
 47     target->transportabilities = new transportabilitiesMap( );
 48     return( 0 );
 49 }
 50 /*
 51 ************************************************************
 52 */
 53 MCGIDI_target_heated *MCGIDI_target_heated_newRead( statusMessageReporting *smr, const char *fileName ) {
 54 
 55     MCGIDI_target_heated *target;
 56 
 57     if( ( target = MCGIDI_target_heated_new( smr ) ) == NULL ) return( NULL );
 58     if( MCGIDI_target_heated_read( smr, target, fileName ) != 0 ) target = (MCGIDI_target_heated *) smr_freeMemory( (void **) &target );
 59     return( target );
 60 }
 61 /*
 62 ************************************************************
 63 */
 64 MCGIDI_target_heated *MCGIDI_target_heated_free( statusMessageReporting *smr, MCGIDI_target_heated *target ) {
 65 
 66     MCGIDI_target_heated_release( smr, target );
 67     smr_freeMemory( (void **) &target );
 68     return( NULL );
 69 }
 70 /*
 71 ************************************************************
 72 */
 73 int MCGIDI_target_heated_release( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
 74 
 75     int ir;
 76 
 77     ptwXY_free( target->crossSection );
 78     ptwX_free( target->crossSectionGrouped );
 79     ptwX_free( target->crossSectionGroupedForSampling );
 80     for( ir = 0; ir < target->numberOfReactions; ir++ ) MCGIDI_reaction_release( NULL, &(target->reactions[ir]) );
 81     smr_freeMemory( (void **) &(target->reactions) );
 82     MCGIDI_POPs_release( &(target->pops) );
 83     smr_freeMemory( (void **) &(target->path) );
 84     smr_freeMemory( (void **) &(target->absPath) );
 85     xDataTOMAL_release( &(target->attributes) );
 86     delete target->transportabilities;
 87     return( 0 );
 88 }
 89 /*
 90 ************************************************************
 91 */
 92 int MCGIDI_target_heated_read( statusMessageReporting *smr, MCGIDI_target_heated *target, const char *fileName ) {
 93 /*
 94 *   If a target has already been read into this target, user must have called MCGIDI_target_heated_release before calling this routine.
 95 *   Otherwise, there will be memory leaks.
 96 */
 97     int n, ir;
 98     xDataTOM_TOM *doc = NULL;
 99     xDataTOM_element *element, *child, *particles, *particleAliases;
100     char const *name, *version, *temperatureStr;
101     char *e1;
102     MCGIDI_reaction *reaction;
103     double crossSectionInit[4] = { 0., 0., 0., 0., };
104     nfu_status status;
105     ptwXYPoints *crossSection;
106 int subtag1_Notice = 0;
107 
108     if( ( target->path = smr_allocateCopyString2( smr, fileName, "path" ) ) == NULL ) goto err;
109     if( ( target->absPath = xDataTOMMisc_getAbsPath( smr, fileName ) ) == NULL ) goto err;
110     if( ( doc = xDataTOM_importFile( smr, fileName ) ) == NULL ) goto err;
111     element = xDataTOM_getDocumentsElement( doc );
112     if( ( version = xDataTOM_getAttributesValueInElement( element, "version" ) ) == NULL ) {
113             smr_setReportError2( smr, smr_unknownID, 1, "version attribute missing from element '%s'", element->name );
114             goto err; }
115     else {
116         if( strcmp( version, "GND 1.3" ) != 0 ) {
117             smr_setReportError2( smr, smr_unknownID, 1, "Unsupported version '%s' for element %s", version, element->name );
118             goto err;
119         }
120     }
121     if( strcmp( element->name, "reactionSuite" ) != 0 ) {
122         smr_setReportError2( smr, smr_unknownID, 1, "input file's top element must be reactionSuite and not %s", element->name );
123         goto err; }
124     else {
125         xDataTOMAL_copyAttributionList( smr, &(target->attributes), &(element->attributes) );
126         particleAliases = xDataTOME_getOneElementByName( smr, element, "aliases", 0 );
127         if( ( particles = xDataTOME_getOneElementByName( smr, element, "particles", 1 ) ) == NULL ) goto err;
128         if( MCGIDI_target_heated_parsePOPs( smr, target, particles, particleAliases ) != 0 ) goto err;
129 
130         if( ( temperatureStr = MCGIDI_misc_pointerToTOMAttributeIfAllOk3( smr, target->absPath, 1, &(target->attributes), "temperature" ) ) == NULL ) goto err;
131         target->temperature_MeV = strtod( temperatureStr, &e1 );
132         while( isspace( *e1 ) ) ++e1; // Loop checking, 11.06.2015, T. Koi
133         target->temperature_MeV *= MCGIDI_misc_getUnitConversionFactor( smr, e1, "MeV/k" );
134         if( !smr_isOk( smr ) ) goto err;
135 
136         if( ( name = MCGIDI_misc_pointerToTOMAttributeIfAllOk3( smr, target->absPath, 1, &(target->attributes), "projectile" ) ) != NULL )
137             target->projectilePOP = MCGIDI_POPs_findParticle( &(target->pops), name );
138         if( !smr_isOk( smr ) ) goto err;
139 
140         if( ( name = MCGIDI_misc_pointerToTOMAttributeIfAllOk3( smr, target->absPath, 1, &(target->attributes), "target" ) ) != NULL )
141         if( !smr_isOk( smr ) ) goto err;
142             target->targetPOP = MCGIDI_POPs_findParticle( &(target->pops), name );
143 
144         n = xDataTOM_numberOfElementsByName( smr, element, "reaction" );
145         if( n == 0 ) {
146             smr_setReportError2( smr, smr_unknownID, 1, "target does not have any reactions: file = '%s'", fileName );
147             goto err;
148         }
149         if( ( target->reactions = (MCGIDI_reaction *) smr_malloc2( smr, n * sizeof( MCGIDI_reaction ), 1, "target->reactions" ) ) == NULL ) goto err;
150 
151         for( ir = 0, child = xDataTOME_getFirstElement( element ); child != NULL; child = xDataTOME_getNextElement( child ) ) {
152             if( strcmp( child->name, "particles" ) == 0 ) continue;
153             if( strcmp( child->name, "styles" ) == 0 ) continue;
154             if( strcmp( child->name, "documentations" ) == 0 ) continue;
155             if( strcmp( child->name, "resonances" ) == 0 ) continue;
156             if( strcmp( child->name, "summedReaction" ) == 0 ) continue;
157             if( strcmp( child->name, "fissionComponent" ) == 0 ) continue;
158             if( strcmp( child->name, "reaction" ) == 0 ) {
159                 double EMin, EMax;
160 
161                 reaction =  &(target->reactions[ir]);
162                 if( MCGIDI_target_heated_parseReaction( smr, child, target, &(target->pops), reaction ) ) goto err;
163                 if( MCGIDI_reaction_getDomain( smr, reaction, &EMin, &EMax ) ) goto err;
164                 if( ir == 0 ) { target->EMin = EMin; target->EMax = EMax; }
165                 if( EMin < target->EMin ) target->EMin = EMin;
166                 if( EMax > target->EMax ) target->EMax = EMax;
167                 for( transportabilitiesMap::const_iterator iter = reaction->transportabilities->begin( ); 
168                             iter != reaction->transportabilities->end( ); ++iter ) {
169                     MCGIDI_misc_updateTransportabilitiesMap( target->transportabilities, iter->first, iter->second );
170                 }
171                 ir++; }
172             else if( strcmp( child->name, "production" ) == 0 ) {
173                 continue; }
174             else if( strcmp( child->name, "aliases" ) == 0 ) {
175                 continue; }
176             else if( strcmp( child->name, "partialGammaProduction" ) == 0 ) {
177                 if( subtag1_Notice == 0 ) printf( "Unsupported reactionSuite sub-tag = '%s'\n", child->name );
178                 subtag1_Notice++; }
179             else {
180                 printf( "Unsupported reactionSuite sub-tag = '%s'\n", child->name );
181             }
182         }
183         crossSectionInit[0] = target->EMin;
184         crossSectionInit[2] = target->EMax;
185         if( ( target->crossSection = ptwXY_create( ptwXY_interpolationLinLin, NULL, 2., 1e-3, 2, 10, 2, crossSectionInit, &status, 0 ) ) == NULL ) {
186             smr_setReportError2( smr, smr_unknownID, 1, "ptwXY_create err = %d: %s\n", status, nfu_statusMessage( status ) );
187             goto err;
188         }
189         for( ir = 0; ir < target->numberOfReactions; ir++ ) {
190             reaction =  &(target->reactions[ir]);
191             if( MCGIDI_reaction_fixDomains( smr, reaction, target->EMin, target->EMax, &status ) ) {
192                 smr_setReportError2( smr, smr_unknownID, 1, "ptwXY_add_ptwXY err = %d: %s\n", status, nfu_statusMessage( status ) );
193                 goto err;
194             }
195             if( ( crossSection = ptwXY_add_ptwXY( target->crossSection, reaction->crossSection, &status ) ) == NULL ) {
196                 smr_setReportError2( smr, smr_unknownID, 1, "ptwXY_add_ptwXY err = %d: %s\n", status, nfu_statusMessage( status ) );
197                 goto err;
198             }
199             target->crossSection = ptwXY_free( target->crossSection );
200             target->crossSection = crossSection;
201         }
202     }
203     xDataTOM_freeTOM( smr, &doc );
204     return( 0 );
205 
206 err:
207     smr_setReportError2( smr, smr_unknownID, 1, "Sub-error while reading file '%s'", fileName );
208     if( doc != NULL ) xDataTOM_freeTOM( smr, &doc );
209     MCGIDI_target_heated_release( smr, target );
210     return( 1 );
211 }
212 /*
213 ************************************************************
214 */
215 static int MCGIDI_target_heated_parsePOPs( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element,
216         xDataTOM_element *particleAliases ) {
217 
218     xDataTOM_element *child;
219 
220     for( child = xDataTOME_getFirstElement( element ); child != NULL; child = xDataTOME_getNextElement( child ) ) {
221         if( strcmp( child->name, "particle" ) ) {
222             smr_setReportError2( smr, smr_unknownID, 1, "invalid element '%s' in %s", child->name, element->name );
223             goto err;
224         }
225         if( MCGIDI_target_heated_parseParticle( smr, target, child, particleAliases ) ) goto err;
226     }
227     return( 0 );
228 
229 err:
230     return( 1 );
231 }
232 /*
233 ************************************************************
234 */
235 static int MCGIDI_target_heated_parseParticle( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element,
236         xDataTOM_element *particleAliases ) {
237 /*
238     This routine, MCGIDI_target_heated_parseParticleLevel, MCGIDI_target_heated_parseParticleGammas handle the parsing of a
239     particle which can have one of the following three forms:
240 
241 1)
242     <particle name="Sc47" genre="nucleus" mass="46.9524440292728 amu"/>
243 2)
244     <particle name="Sc46" genre="nucleus" mass="45.9551770270807 amu">
245         <level name="Sc46_e0" label="0" energy="0 eV"/>
246         <level name="Sc46_e1" label="1" energy="142500 eV"/></particle>
247 3)
248     <particle name="S36" genre="nucleus" mass="35.9669735654569 amu">
249         <level name="S36_e0" label="0" energy="0 eV" spin="0"/>
250         <level name="S36_e1" label="1" energy="3.291e6 eV">
251             <gamma finalLevel="S36_e0" probability="1.0"/></level>
252         <level name="S36_e2" label="2" energy="3.346e6 eV">
253             <gamma finalLevel="S36_e0" probability="1.0"/></level>
254         <level name="S36_e3" label="3" energy="4191990 eV">
255             <gamma finalLevel="S36_e1" probability="1.0"/></level>
256         <level name="S36_e4" label="4" energy="4522990 eV">
257             <gamma finalLevel="S36_e1" probability="0.248"/>
258             <gamma finalLevel="S36_e0" probability="0.752"/></level>
259         <level name="S36_e5" label="5" energy="4574990 eV">
260             <gamma finalLevel="S36_e1" probability="1.0"/></level>
261         <level name="S36_c" label="c" energy="u:4999990 eV"/></particle>
262 */
263     int globalParticle = 1;
264     char const *name = NULL, *mass = NULL;      /* Do not free name or mass, do not own them. */
265     double mass_MeV;
266     xDataTOM_element *child;
267     MCGIDI_POP *pop;
268 
269     if( ( name = xDataTOM_getAttributesValueInElement( element, "name" ) ) == NULL ) {
270         smr_setReportError2p( smr, smr_unknownID, 1, "particle missing name attribute" );
271         goto err;
272     }
273     if( ( mass = xDataTOM_getAttributesValueInElement( element, "mass" ) ) == NULL ) {
274         smr_setReportError2( smr, smr_unknownID, 1, "particle '%s' missing mass attribute", name );
275         goto err;
276     }
277     if( MCGIDI_misc_PQUStringToDouble( smr, mass, "amu", MCGIDI_AMU2MeV, &mass_MeV ) ) goto err;
278     if( ( pop = MCGIDI_POPs_addParticleIfNeeded( smr, &(target->pops), name, mass_MeV, 0., NULL, globalParticle ) ) == NULL ) goto err;
279 
280     for( child = xDataTOME_getFirstElement( element ); child != NULL; child = xDataTOME_getNextElement( child ) ) {
281         if( strcmp( child->name, "level" ) ) {
282             smr_setReportError2( smr, smr_unknownID, 1, "invalid element '%s' in %s", child->name, element->name );
283             goto err;
284         }
285         if( MCGIDI_target_heated_parseParticleLevel( smr, target, child, pop, mass_MeV, particleAliases ) ) goto err;
286     }
287 
288     return( 0 );
289 
290 err:
291     return( 1 );
292 }
293 /*
294 ************************************************************
295 */
296 static int MCGIDI_target_heated_parseParticleLevel( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element, MCGIDI_POP *parent,
297     double mass_MeV, xDataTOM_element *particleAliases ) {
298 
299     int globalParticle = 0;
300     char const *name, *level, *aliasValue;     /* Do not free any of these as they are owned by called routine. */
301     double level_MeV = 0.;
302     xDataTOM_element *alias;
303 
304     if( ( name = xDataTOM_getAttributesValueInElement( element, "name" ) ) == NULL ) {
305         smr_setReportError2p( smr, smr_unknownID, 1, "particle missing name attribute" );
306         goto err;
307     }
308     if( ( level = xDataTOM_getAttributesValueInElement( element, "energy" ) ) == NULL ) {
309         smr_setReportError2( smr, smr_unknownID, 1, "particle '%s' level missing energy attribute", name );
310         goto err;
311     }
312         /* Special case for 'c' labels. Correct mass is only needed for two-body. */
313     if( level[0] != 'u' ) if( MCGIDI_misc_PQUStringToDoubleInUnitOf( smr, level, "MeV", &level_MeV ) ) goto err;
314     for( alias = xDataTOME_getFirstElement( particleAliases ); alias != NULL; alias = xDataTOME_getNextElement( alias ) ) {
315         if( ( aliasValue = xDataTOM_getAttributesValueInElement( alias, "value" ) ) == NULL ) {
316             smr_setReportError2p( smr, smr_unknownID, 1, "particle missing name attribute" );
317             goto err;
318         }
319         if( strcmp( aliasValue, name ) == 0 ) globalParticle = 1;
320     }
321     if( MCGIDI_POPs_addParticleIfNeeded( smr, &(target->pops), name, mass_MeV + level_MeV, level_MeV, parent, globalParticle ) == NULL ) goto err;
322 
323     return( MCGIDI_target_heated_parseParticleGammas( smr, target, element, name ) );
324 
325 err:
326     return( 1 );
327 }
328 /*
329 ************************************************************
330 */
331 static int MCGIDI_target_heated_parseParticleGammas( statusMessageReporting *smr, MCGIDI_target_heated *target, xDataTOM_element *element, char const *name ) {
332 
333     int gammaCounts = 0;
334     MCGIDI_POP *pop = MCGIDI_POPs_findParticle( &(target->pops), name );
335     xDataTOM_element *child;
336     MCGIDI_GammaBranching *gammas = NULL;
337     char const *finalLevelString;
338     double probability;
339 
340     for( child = xDataTOME_getFirstElement( element ); child != NULL; child = xDataTOME_getNextElement( child ) ) {
341         if( strcmp( child->name, "gamma" ) ) {
342             smr_setReportError2( smr, smr_unknownID, 1, "invalid element '%s' in %s", child->name, element->name );
343             goto err;
344         }
345         gammaCounts++;
346     }
347     if( gammaCounts > 0 ) {
348         if( ( gammas = (MCGIDI_GammaBranching *) smr_malloc2( smr, gammaCounts * sizeof( MCGIDI_GammaBranching), 0, "gammas" ) ) == NULL ) goto err;
349         for( child = xDataTOME_getFirstElement( element ); child != NULL; child = xDataTOME_getNextElement( child ) ) {
350             if( ( finalLevelString = xDataTOM_getAttributesValueInElement( child, "finalLevel" ) ) == NULL ) {
351                 smr_setReportError2p( smr, smr_unknownID, 1, "gamma missing 'finalLevel'" );
352                 goto err;
353             }
354             if( xDataTOME_convertAttributeToDouble( smr, child, "probability", &probability ) != 0 ) {
355                 smr_setReportError2p( smr, smr_unknownID, 1, "gamma missing 'probability' attribute" );
356                 goto err;
357             }
358         }
359     }
360     pop->numberOfGammaBranchs = gammaCounts;
361     pop->gammas = gammas;
362 
363     return( 0 );
364 
365 err:
366     if( gammas != NULL ) smr_freeMemory( (void **) &gammas );
367     return( 1 );
368 }
369 /*
370 ************************************************************
371 */
372 static int MCGIDI_target_heated_parseReaction( statusMessageReporting *smr, xDataTOM_element *element, MCGIDI_target_heated *target, 
373         MCGIDI_POPs *pops, MCGIDI_reaction *reaction ) {
374 
375     if( MCGIDI_reaction_parseFromTOM( smr, element, target, pops, reaction ) ) goto err;
376     target->numberOfReactions++;
377     
378     return( 0 );
379 
380 err:
381     smr_setReportError2( smr, smr_unknownID, 1, "%s\n", xDataTOM_getAttributesValueInElement( element, "outputChannel" ) );
382     return( 1 );
383 }
384 /*
385 ************************************************************
386 */
387 int MCGIDI_target_heated_numberOfReactions( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
388 
389     return( target->numberOfReactions );
390 }
391 /*
392 ************************************************************
393 */
394 int MCGIDI_target_heated_numberOfProductionReactions( statusMessageReporting * /*smr*/, MCGIDI_target_heated * /*target*/ ) {
395 
396     return( 0 );
397 }
398 /*
399 ************************************************************
400 */
401 MCGIDI_reaction *MCGIDI_target_heated_getReactionAtIndex( MCGIDI_target_heated *target, int index ) {
402 
403     if( ( index >= 0 ) && ( index < target->numberOfReactions ) ) return( &(target->reactions[index]) );
404     return( NULL );
405 }
406 /*
407 ************************************************************
408 */
409 MCGIDI_reaction *MCGIDI_target_heated_getReactionAtIndex_smr( statusMessageReporting *smr, MCGIDI_target_heated *target, int index ) {
410 
411     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex( target, index );
412 
413     if( reaction == NULL ) {
414         smr_setReportError2( smr, smr_unknownID, 1, "bad reaction index = %d for %s + %s", index, target->projectilePOP->name, target->targetPOP->name ); 
415     }
416     return( reaction );
417 }
418 #if 0
419 /*
420 ************************************************************
421 */
422 MCGIDI_channel *MCGIDI_target_heated_getProductionReactionAtIndex( MCGIDI_target_heated *target, int index ) {
423 
424     MCGIDI_channel *channel = NULL;
425 
426     if( ( index >= 0 ) && ( index < target->nProductionReactions ) ) channel = target->productionReactions[index];
427     return( channel );
428 }
429 #endif
430 /*
431 ************************************************************
432 */
433 MCGIDI_POP *MCGIDI_target_heated_getPOPForProjectile( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
434 
435     return( target->projectilePOP );
436 }
437 /*
438 ************************************************************
439 */
440 MCGIDI_POP *MCGIDI_target_heated_getPOPForTarget( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
441 
442     return( target->targetPOP );
443 }
444 /*
445 ************************************************************
446 */
447 double MCGIDI_target_heated_getProjectileMass_MeV( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
448 
449     return( MCGIDI_POP_getMass_MeV( target->projectilePOP ) );
450 }
451 /*
452 ************************************************************
453 */
454 double MCGIDI_target_heated_getTargetMass_MeV( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
455 
456     return( MCGIDI_POP_getMass_MeV( target->targetPOP ) );
457 }
458 /*
459 ************************************************************
460 */
461 double MCGIDI_target_heated_getTotalCrossSectionAtE( statusMessageReporting *smr, MCGIDI_target_heated *target, 
462         MCGIDI_quantitiesLookupModes &modes, bool sampling ) {
463 
464     double xsec;
465 
466     if( modes.getCrossSectionMode( ) == MCGIDI_quantityLookupMode_pointwise ) {
467         double e_in = modes.getProjectileEnergy( );
468 
469         if( e_in < target->EMin ) e_in = target->EMin;
470         if( e_in > target->EMax ) e_in = target->EMax;
471         ptwXY_getValueAtX( target->crossSection, e_in, &xsec ); }
472     else if( modes.getCrossSectionMode( ) == MCGIDI_quantityLookupMode_grouped ) {
473         int index = modes.getGroupIndex( );
474         double *xSecP;
475 
476         if( sampling ) {
477             xSecP = ptwX_getPointAtIndex( target->crossSectionGroupedForSampling, index ); }
478         else {
479             xSecP = ptwX_getPointAtIndex( target->crossSectionGrouped, index );
480         }
481 
482         if( xSecP != NULL ) {
483             xsec = *xSecP; }
484         else {
485             xsec = 0.;
486             smr_setReportError2( smr, smr_unknownID, 1, "Invalid cross section group index %d", index, (int) ptwX_length( target->crossSectionGrouped ) );
487         } }
488     else {
489         xsec = 0.;
490     }
491     return( xsec );
492 }
493 /*
494 ************************************************************
495 */
496 double MCGIDI_target_heated_getIndexReactionCrossSectionAtE( statusMessageReporting *smr, MCGIDI_target_heated *target, int index, 
497         MCGIDI_quantitiesLookupModes &modes, bool sampling ) {
498 
499     double xsec = 0.;
500     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex_smr( smr, target, index );
501 
502     if( reaction != NULL ) xsec = MCGIDI_reaction_getCrossSectionAtE( smr, reaction, modes, sampling );
503     return( xsec );
504 }
505 /*
506 ************************************************************
507 */
508 int MCGIDI_target_heated_sampleIndexReactionProductsAtE( statusMessageReporting *smr, MCGIDI_target_heated *target, int index, 
509         MCGIDI_quantitiesLookupModes &modes, MCGIDI_decaySamplingInfo *decaySamplingInfo, MCGIDI_sampledProductsDatas *productDatas ) {
510 
511     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex_smr( smr, target, index );
512 
513     productDatas->numberOfProducts = 0;
514     if( reaction == NULL ) return( -1 );
515     return( MCGIDI_outputChannel_sampleProductsAtE( smr, &(reaction->outputChannel), modes, decaySamplingInfo, productDatas, NULL ) );
516 }
517 /*
518 ************************************************************
519 */
520 double MCGIDI_target_heated_getReactionsThreshold( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target, int index ) {
521 
522     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex( target, index );
523 
524     if( reaction == NULL ) return( -1 );
525     return( reaction->EMin );
526 }
527 /*
528 ************************************************************
529 */
530 int MCGIDI_target_heated_getReactionsDomain( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target, int index, double *EMin, double *EMax ) {
531 
532     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex( target, index );
533 
534     if( reaction == NULL ) return( -1 );
535     *EMin = reaction->EMin;
536     *EMax = reaction->EMax;
537     return( 0 );
538 }
539 /*
540 ************************************************************
541 */
542 double MCGIDI_target_heated_getIndexReactionFinalQ( statusMessageReporting *smr, MCGIDI_target_heated *target, int index, 
543         MCGIDI_quantitiesLookupModes &modes ) {
544 
545     MCGIDI_reaction *reaction = MCGIDI_target_heated_getReactionAtIndex_smr( smr, target, index );
546 
547     if( reaction == NULL ) return( 0. );
548     return( MCGIDI_reaction_getFinalQ( smr, reaction, modes ) );
549 }
550 /*
551 ************************************************************
552 */
553 std::map<int, enum MCGIDI_transportability> const *MCGIDI_target_heated_getUniqueProducts( statusMessageReporting * /*smr*/, MCGIDI_target_heated *target ) {
554 
555     return( target->transportabilities );
556 }
557 /*
558 ************************************************************
559 */
560 int MCGIDI_target_heated_recast( statusMessageReporting *smr, MCGIDI_target_heated *target, GIDI_settings &settings ) {
561 
562     int ir, projectilePoPID = target->projectilePOP->globalPoPsIndex;
563     ptwXPoints *totalGroupedCrossSection = NULL;
564     GIDI_settings_particle const *projectileSettings = settings.getParticle( projectilePoPID );
565     nfu_status status_nf;
566 
567     if( projectileSettings == NULL ) {
568         smr_setReportError2( smr, smr_unknownID, 1, "Settings missing for projectile %s", target->projectilePOP->name ); 
569         return( 1 );
570     }
571     target->crossSectionGrouped = ptwX_free( target->crossSectionGrouped );
572     target->crossSectionGroupedForSampling = ptwX_free( target->crossSectionGroupedForSampling );
573     if( projectileSettings->isEnergyMode_grouped( ) ) {
574         int64_t numberOfGroups = projectileSettings->getNumberOfGroups( );
575 
576         if( ( totalGroupedCrossSection = ptwX_createLine( numberOfGroups, numberOfGroups, 0, 0, &status_nf ) ) == NULL ) {
577             smr_setReportError2( smr, smr_unknownID, 1, "totalGroupedCrossSection allocation failed: status_nf = %d, '%s'", 
578                     status_nf, nfu_statusMessage( status_nf ) ); 
579             goto err;
580         }
581     }
582 
583     for( ir = 0; ir < target->numberOfReactions; ++ir ) {
584         if( MCGIDI_reaction_recast( smr, &(target->reactions[ir]), settings, projectileSettings, target->temperature_MeV, totalGroupedCrossSection ) != 0 ) goto err;
585     }
586     if( projectileSettings->isEnergyMode_grouped( ) ) {
587         if( ( target->crossSectionGroupedForSampling = ptwX_clone( totalGroupedCrossSection, &status_nf ) ) == NULL ) {
588             smr_setReportError2( smr, smr_unknownID, 1, "totalGroupedCrossSection allocation failed: status_nf = %d, '%s'", 
589                     status_nf, nfu_statusMessage( status_nf ) ); 
590             goto err;
591         }
592         for( ir = 0; ir < target->numberOfReactions; ++ir ) {
593             int index = target->reactions[ir].thresholdGroupIndex;
594 
595             if( index > -1 ) {
596                 double xSec = target->reactions[ir].thresholdGroupedDeltaCrossSection +
597                         ptwX_getPointAtIndex_Unsafely( target->crossSectionGroupedForSampling, index );
598 
599                 ptwX_setPointAtIndex( target->crossSectionGroupedForSampling, index, xSec );
600             }
601         }
602     }
603     target->crossSectionGrouped = totalGroupedCrossSection;
604     totalGroupedCrossSection = NULL;
605 
606     return( 0 );
607 
608 err:
609     ptwX_free( totalGroupedCrossSection );
610     target->crossSectionGroupedForSampling = ptwX_free( target->crossSectionGroupedForSampling );
611     return( 1 );
612 }
613 
614 #if defined __cplusplus
615 }
616 #endif
617