Geant4  10.00.p02
tpia_angularEnergy.cc
Go to the documentation of this file.
1 /*
2 # <<BEGIN-copyright>>
3 # Copyright (c) 2010, Lawrence Livermore National Security, LLC.
4 # Produced at the Lawrence Livermore National Laboratory
5 # Written by Bret R. Beck, beck6@llnl.gov.
6 # CODE-461393
7 # All rights reserved.
8 #
9 # This file is part of GIDI. For details, see nuclear.llnl.gov.
10 # Please also read the "Additional BSD Notice" at nuclear.llnl.gov.
11 #
12 # Redistribution and use in source and binary forms, with or without modification,
13 # are permitted provided that the following conditions are met:
14 #
15 # 1) Redistributions of source code must retain the above copyright notice,
16 # this list of conditions and the disclaimer below.
17 # 2) Redistributions in binary form must reproduce the above copyright notice,
18 # this list of conditions and the disclaimer (as noted below) in the
19 # documentation and/or other materials provided with the distribution.
20 # 3) Neither the name of the LLNS/LLNL nor the names of its contributors may be
21 # used to endorse or promote products derived from this software without
22 # specific prior written permission.
23 #
24 # THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
25 # EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
26 # OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
27 # SHALL LAWRENCE LIVERMORE NATIONAL SECURITY, LLC, THE U.S. DEPARTMENT OF ENERGY OR
28 # CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
29 # CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
30 # OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
31 # AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
32 # (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE,
33 # EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
34 # <<END-copyright>>
35 */
36 #include <string.h>
37 
38 #include <tpia_target.h>
39 #include <tpia_misc.h>
40 
41 #if defined __cplusplus
42 namespace GIDI {
43 using namespace GIDI;
44 #endif
45 
46 /*
47 ************************************************************
48 */
49 int tpia_angularEnergy_initialize( statusMessageReporting *smr, tpia_angularEnergy *angularEnergy ) {
50 
51  memset( angularEnergy, 0, sizeof( tpia_angularEnergy ) );
52  if( tpia_frame_setFromString( smr, "", "", 0, &(angularEnergy->frame) ) ) return( 1 );
53  angularEnergy->binned.numberOfEs = 0;
54  angularEnergy->binned.energies = NULL;
55  return( 0 );
56 }
57 /*
58 ************************************************************
59 */
60 int tpia_angularEnergy_release( statusMessageReporting *smr, tpia_angularEnergy *angularEnergy ) {
61 
62  int i;
63 
64  for( i = 0; i < angularEnergy->binned.numberOfEs; i++ ) xData_free( smr, angularEnergy->binned.energies[i].energies );
65  //angularEnergy->binned.energies = xData_free( smr, angularEnergy->binned.energies );
66  angularEnergy->binned.energies = (tpia_EqualProbableBinSpectra*) xData_free( smr, angularEnergy->binned.energies );
67  tpia_angularEnergy_initialize( smr, angularEnergy );
68  return( 0 );
69 }
70 /*
71 ************************************************************
72 */
73 int tpia_angularEnergy_getFromElement( statusMessageReporting *smr, xData_element *angularEnergyElement, tpia_angularEnergy *angularEnergy ) {
74 
75  int i, j, status = 1;
76  xData_Int nBins, size, n, index;
77  xData_element *epbElement, *element;
78  xData_elementList *list;
79  tpia_EqualProbableBinSpectra *mus;
80 
81  xData_addToAccessed( smr, angularEnergyElement, 1 );
82  if( ( tpia_frame_setFromElement( smr, angularEnergyElement, 4, &(angularEnergy->frame) ) ) != 0 ) return( 1 );
83  //if( ( epbElement = xData_getOneElementByTagName( smr, angularEnergyElement, "equalProbableBins", 0 ) ) == NULL ) return( 1 );
84  if( ( epbElement = xData_getOneElementByTagName( smr, angularEnergyElement, (char*)"equalProbableBins", 0 ) ) == NULL ) return( 1 );
85  xData_addToAccessed( smr, epbElement, 1 );
86  if( xData_convertAttributeTo_xData_Int( smr, epbElement, "nBins", &nBins ) != 0 ) {
87  tpia_misc_setMessageError_Element( smr, NULL, epbElement, __FILE__, __LINE__, 1, "missing or invalid nBins attribute" );
88  return( 1 );
89  }
90  list = xData_getElementsByTagNameAndSort( smr, epbElement, "energy", NULL, NULL );
91  if( list->n == 0 ) {
92  tpia_misc_setMessageError_Element( smr, NULL, epbElement, __FILE__, __LINE__, 1, "bins does not contain any energy elements" ); }
93  else {
94  size = list->n * sizeof( tpia_EqualProbableBinSpectra );
95  angularEnergy->binned.nBins = (int) nBins;
96  //if( ( angularEnergy->binned.energies = xData_malloc2( smr, size, 1, "energies" ) ) != NULL ) {
97  if( ( angularEnergy->binned.energies = (tpia_EqualProbableBinSpectra*) xData_malloc2( smr, size, 1, "energies" ) ) != NULL ) {
98  status = 0;
99  for( i = 0, mus = angularEnergy->binned.energies; i < list->n; i++, mus++ ) {
100  mus->iValue = 0;
101  element = list->items[i].element;
102  if( xData_convertAttributeTo_xData_Int( smr, element, "index", &index ) != 0 ) {
103  tpia_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "missing or invalid index attribute" );
104  status = 1;
105  break;
106  }
107  if( index != i ) {
108  tpia_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "index = %lld is not incremental", index );
109  status = 1;
110  break;
111  }
112  if( ( j = xData_convertAttributeToDouble( smr, element, "value", &(mus->dValue) ) ) != 0 ) {
113  if( j == 1 ) {
114  tpia_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "element does not have value attribute" ); }
115  else {
116  tpia_misc_setMessageError_Element( smr, NULL, element, __FILE__, __LINE__, 1, "failed to convert value attribute to double" );
117  }
118  status = 1;
119  break;
120  }
121  if( ( mus->energies = tpia_misc_getEqualProbableBins( smr, element, "mu", nBins, &n ) ) == NULL ) {
122  status = 1;
123  break; }
124  else {
125  mus->numberOfEs = n;
126  }
127  angularEnergy->binned.numberOfEs++;
128  }
129  }
130  }
131  xData_freeElementList( smr, list );
132  return( status );
133 }
134 /*
135 ************************************************************
136 */
137 int tpia_angularEnergy_SampleEp( statusMessageReporting *smr, tpia_angularEnergy *angularEnergy, tpia_decaySamplingInfo *decaySamplingInfo ) {
138 
139  int iE1, iE2;
140  tpia_EqualProbableBinSpectra *energies = angularEnergy->binned.energies;
141  double Ep, Ep1, Ep2, f, e_in = decaySamplingInfo->e_in;
142 /*
143 Currently, only equal probable binning is supported.
144 Need to return frame info for Ep also.
145 */
146 
147  if( !smr_isOk( smr ) ) return( 1 );
148  if( angularEnergy->binned.numberOfEs == 0 ) return( 1 );
149  for( iE2 = 0; iE2 < angularEnergy->binned.numberOfEs; iE2++ ) if( energies[iE2].dValue >= e_in ) break;
150  if( iE2 == 0 ) {
151  iE1 = iE2; }
152  else if( iE2 == angularEnergy->binned.numberOfEs ) {
153  iE1 = iE2 = angularEnergy->binned.numberOfEs - 1; }
154  else {
155  iE1 = iE2 - 1;
156  }
157  tpia_misc_sampleEqualProbableBin( smr, decaySamplingInfo, decaySamplingInfo->mu, angularEnergy->binned.nBins, &(energies[iE1]), &Ep1 );
158  if( iE1 == iE2 ) {
159  Ep = Ep1; }
160  else {
161  tpia_misc_sampleEqualProbableBin( smr, decaySamplingInfo, decaySamplingInfo->mu, angularEnergy->binned.nBins, &(energies[iE2]), &Ep2 );
162  f = ( energies[iE2].dValue - e_in ) / ( energies[iE2].dValue - energies[iE1].dValue );
163  Ep = f * Ep1 + ( 1 - f ) * Ep2;
164  }
165  decaySamplingInfo->Ep = Ep;
166 
167  return( 0 );
168 }
169 
170 #if defined __cplusplus
171 }
172 #endif
xData_Int xData_convertAttributeTo_xData_Int(statusMessageReporting *smr, xData_element *element, const char *name, xData_Int *n)
Definition: xData.cc:710
xData_element * xData_getOneElementByTagName(statusMessageReporting *smr, xData_element *element, char *name, int required)
Definition: xData.cc:810
int xData_addToAccessed(statusMessageReporting *, xData_element *element, int increment)
Definition: xData.cc:1078
int tpia_misc_setMessageError_Element(statusMessageReporting *smr, void *userInterface, xData_element *element, const char *file, int line, int code, const char *fmt,...)
Definition: tpia_misc.cc:183
int tpia_angularEnergy_initialize(statusMessageReporting *smr, tpia_angularEnergy *angularEnergy)
int tpia_angularEnergy_release(statusMessageReporting *smr, tpia_angularEnergy *angularEnergy)
int tpia_frame_setFromString(statusMessageReporting *smr, const char *forItem, const char *value, int dimension, tpia_data_frame *frame)
Definition: tpia_frame.cc:75
xData_elementList * xData_getElementsByTagNameAndSort(statusMessageReporting *smr, xData_element *element, const char *tagName, const char *sortAttributeName, xData_sortElementFunc sortFunction)
Definition: xData.cc:788
int smr_isOk(statusMessageReporting *smr)
tpia_EqualProbableBinSpectrum * tpia_misc_getEqualProbableBins(statusMessageReporting *smr, xData_element *parent, const char *name, xData_Int nBins, xData_Int *n)
Definition: tpia_misc.cc:333
int tpia_angularEnergy_SampleEp(statusMessageReporting *smr, tpia_angularEnergy *angularEnergy, tpia_decaySamplingInfo *decaySamplingInfo)
const G4int n
int tpia_frame_setFromElement(statusMessageReporting *smr, xData_element *element, int dimension, tpia_data_frame *frame)
Definition: tpia_frame.cc:62
int tpia_angularEnergy_getFromElement(statusMessageReporting *smr, xData_element *angularEnergyElement, tpia_angularEnergy *angularEnergy)
void * xData_free(statusMessageReporting *, void *p)
Definition: xDataMisc.cc:89
int tpia_misc_sampleEqualProbableBin(statusMessageReporting *, tpia_decaySamplingInfo *decaySamplingInfo, double e_in, int nBins, tpia_EqualProbableBinSpectra *binned, double *value_)
Definition: tpia_misc.cc:420
int xData_convertAttributeToDouble(statusMessageReporting *smr, xData_element *element, const char *name, double *d)
Definition: xData.cc:730
void xData_freeElementList(statusMessageReporting *smr, xData_elementList *list)
Definition: xData.cc:832