| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // neutron_hp -- source file 26 // neutron_hp -- source file
27 // J.P. Wellisch, Nov-1996 27 // J.P. Wellisch, Nov-1996
28 // A prototype of the low energy neutron trans 28 // A prototype of the low energy neutron transport model.
29 // 29 //
30 // 13-Jan-06 fix in Sample by T. Koi 30 // 13-Jan-06 fix in Sample by T. Koi
31 // 31 //
32 // P. Arce, June-2014 Conversion neutron_hp to 32 // P. Arce, June-2014 Conversion neutron_hp to particle_hp
33 // 33 //
34 #include "G4ParticleHPPartial.hh" << 34 #include "G4ParticleHPPartial.hh"
35 <<
36 #include "G4ParticleHPInterpolator.hh" 35 #include "G4ParticleHPInterpolator.hh"
37 #include "Randomize.hh" 36 #include "Randomize.hh"
38 37
39 G4ParticleHPVector* G4ParticleHPPartial::GetY( << 38 G4ParticleHPVector * G4ParticleHPPartial::GetY(G4double e1)
40 { << 39 {
41 auto aBuffer = new G4ParticleHPVector(); << 40 G4ParticleHPVector * aBuffer = new G4ParticleHPVector();
42 G4int i; << 41 G4int i;
43 if (nData == 1) { << 42 if(nData==1)
44 for (i = 0; i < data[0].GetVectorLength(); << 43 {
45 aBuffer->SetInterpolationManager(data[0] << 44 for(i=0; i<data[0].GetVectorLength(); i++)
46 aBuffer->SetData(i, data[0].GetX(i), dat << 45 {
>> 46 aBuffer->SetInterpolationManager(data[0].GetInterpolationManager());
>> 47 aBuffer->SetData(i , data[0].GetX(i), data[0].GetY(i));
>> 48 }
>> 49 return aBuffer;
47 } 50 }
48 return aBuffer; << 51 for (i=0; i<nData; i++)
49 } << 52 {
50 for (i = 0; i < nData; i++) { << 53 if(X[i]>e1) break;
51 if (X[i] > e1) break; <<
52 } <<
53 if (i == nData) i--; <<
54 if (0 == i) i = 1; <<
55 G4double x1, x2, y1, y2, y; <<
56 G4int i1 = 0, ib = 0; <<
57 G4double E1 = X[i - 1]; <<
58 G4double E2 = X[i]; <<
59 for (G4int ii = 0; ii < data[i].GetVectorLen <<
60 x1 = data[i - 1].GetX(std::min(i1, data[i <<
61 x2 = data[i].GetX(ii); <<
62 if (x1 < x2 && i1 < data[i - 1].GetVectorL <<
63 y1 = data[i - 1].GetY(i1); <<
64 y2 = data[i].GetY(x1); <<
65 if (E2 - E1 != 0) { <<
66 y = theInt.Interpolate(theManager.GetS <<
67 } <<
68 else { <<
69 y = 0.5 * (y1 + y2); <<
70 } <<
71 aBuffer->SetData(ib, x1, y); <<
72 aBuffer->SetScheme(ib++, data[i - 1].Get <<
73 i1++; <<
74 if (x2 - x1 > 0.001 * x1) { <<
75 ii--; <<
76 } <<
77 } <<
78 else { <<
79 y1 = data[i - 1].GetY(x2); <<
80 y2 = data[i].GetY(ii); <<
81 if (E2 - E1 != 0) { <<
82 y = theInt.Interpolate(theManager.GetS <<
83 } <<
84 else { <<
85 y = 0.5 * (y1 + y2); <<
86 } <<
87 aBuffer->SetData(ib, x2, y); <<
88 aBuffer->SetScheme(ib++, data[i].GetSche <<
89 if (x1 - x2 < 0.001 * x2) i1++; <<
90 } 54 }
91 } << 55 if(i==nData) i--;
92 return aBuffer; << 56 if(0==i) i=1;
93 } << 57 G4double x1,x2,y1,y2,y;
94 << 58 G4int i1=0, ib=0;
95 G4double G4ParticleHPPartial::Sample(G4double << 59 G4double E1 = X[i-1];
96 { << 60 G4double E2 = X[i];
97 G4int i; << 61 for(G4int ii=0; ii<data[i].GetVectorLength(); ii++)
98 for (i = 0; i < nData; i++) { << 62 {
99 if (x < X[i]) break; << 63 x1 = data[i-1].GetX(std::min(i1, data[i-1].GetVectorLength()-1));
100 } << 64 x2 = data[i].GetX(ii);
101 G4ParticleHPVector theBuff; << 65 if(x1<x2&&i1<data[i-1].GetVectorLength())
102 if (i == 0) { << 66 {
103 theBuff.SetInterpolationManager(data[0].Ge << 67 y1 = data[i-1].GetY(i1);
104 for (G4int ii = 0; ii < GetNEntries(0); ii << 68 y2 = data[i].GetY(x1);
105 theBuff.SetX(ii, GetX(0, ii)); << 69 if(E2-E1!=0)
106 theBuff.SetY(ii, GetY(0, ii)); << 70 {
>> 71 y = theInt.Interpolate(theManager.GetScheme(i), e1, E1, E2, y1, y2);
>> 72 }
>> 73 else
>> 74 {
>> 75 y = 0.5*(y1+y2);
>> 76 }
>> 77 aBuffer->SetData(ib, x1, y);
>> 78 aBuffer->SetScheme(ib++, data[i-1].GetScheme(i1));
>> 79 i1++;
>> 80 if(x2-x1>0.001*x1)
>> 81 {
>> 82 ii--;
>> 83 }
>> 84 }
>> 85 else
>> 86 {
>> 87 y1 = data[i-1].GetY(x2);
>> 88 y2 = data[i].GetY(ii);
>> 89 if(E2-E1!=0)
>> 90 {
>> 91 y = theInt.Interpolate(theManager.GetScheme(i), e1, E1, E2, y1, y2);
>> 92 }
>> 93 else
>> 94 {
>> 95 y = 0.5*(y1+y2);
>> 96 }
>> 97 aBuffer->SetData(ib, x2, y);
>> 98 aBuffer->SetScheme(ib++, data[i].GetScheme(ii));
>> 99 if(x1-x2<0.001*x2) i1++;
>> 100 }
107 } 101 }
>> 102 return aBuffer;
108 } 103 }
109 // else if(i==nData-1) this line will be del << 104
110 else if (i == nData) { << 105 G4double G4ParticleHPPartial::Sample(G4double x)
111 for (i = 0; i < GetNEntries(nData - 1); i+ << 106 {
112 theBuff.SetX(i, GetX(nData - 1, i)); << 107 G4int i;
113 theBuff.SetY(i, GetY(nData - 1, i)); << 108 for (i=0; i<nData; i++)
114 theBuff.SetInterpolationManager(data[nDa << 109 {
>> 110 if(x<X[i]) break;
115 } 111 }
116 } << 112 G4ParticleHPVector theBuff;
117 else { << 113 if(i==0)
118 G4int low = i - 1; << 114 {
119 G4int high = low + 1; << 115 theBuff.SetInterpolationManager(data[0].GetInterpolationManager());
120 G4double x1, x2, y1, y2; << 116 for(G4int ii=0;ii<GetNEntries(0);ii++)
121 G4int i1 = 0, i2 = 0, ii = 0; << 117 {
122 x1 = X[low]; << 118 theBuff.SetX(ii, GetX(0,ii));
123 x2 = X[high]; << 119 theBuff.SetY(ii, GetY(0,ii));
124 while (i1 < GetNEntries(low) || i2 < GetNE <<
125 { <<
126 if ((GetX(low, i1) < GetX(high, i2) && i <<
127 theBuff.SetX(ii, GetX(low, i1)); <<
128 y1 = GetY(low, i1); <<
129 y2 = GetY(high, GetX(low, i1)); // pr <<
130 theBuff.SetY(ii, theInt.Interpolate(th <<
131 y2 <<
132 theBuff.SetScheme(ii, data[low].GetSch <<
133 if (std::abs(GetX(low, i1) - GetX(high <<
134 i1++; <<
135 ii++; <<
136 } 120 }
137 else { << 121 }
138 theBuff.SetX(ii, GetX(high, i2)); << 122 //else if(i==nData-1) this line will be delete
139 //************************************ << 123 else if ( i == nData )
140 // Change by E.Mendoza and D.Cano (CIE << 124 {
141 // y1 = GetY(high,i2); << 125 for(i=0;i<GetNEntries(nData-1);i++)
142 // y2 = GetY(low, GetX(high,i2)); //p << 126 {
143 y2 = GetY(high, i2); << 127 theBuff.SetX(i, GetX(nData-1,i));
144 y1 = GetY(low, GetX(high, i2)); // pr << 128 theBuff.SetY(i, GetY(nData-1,i));
145 //************************************ << 129 theBuff.SetInterpolationManager(data[nData-1].GetInterpolationManager());
146 theBuff.SetY(ii, theInt.Interpolate(th << 130 }
147 y2 << 131 }
148 theBuff.SetScheme(ii, data[high].GetSc << 132 else
149 if (std::abs(GetX(low, i1) - GetX(high << 133 {
150 i2++; << 134 G4int low = i-1;
151 ii++; << 135 G4int high = low+1;
>> 136 G4double x1,x2,y1,y2;
>> 137 G4int i1=0, i2=0, ii=0;
>> 138 x1 = X[low];
>> 139 x2 = X[high];
>> 140 while(i1<GetNEntries(low)||i2<GetNEntries(high)) // Loop checking, 11.05.2015, T. Koi
>> 141 {
>> 142 if( (GetX(low,i1)<GetX(high,i2) && i1<GetNEntries(low))
>> 143 ||(i2==GetNEntries(high)) )
>> 144 {
>> 145 theBuff.SetX(ii, GetX(low,i1));
>> 146 y1 = GetY(low,i1);
>> 147 y2 = GetY(high, GetX(low,i1)); //prob at ident theta
>> 148 theBuff.SetY(ii, theInt.Interpolate(theManager.GetScheme(high),
>> 149 x, x1, x2, y1, y2)); //energy interpol
>> 150 theBuff.SetScheme(ii, data[low].GetScheme(i1));
>> 151 if(std::abs(GetX(low,i1)-GetX(high,i2))<0.001) i2++;
>> 152 i1++;
>> 153 ii++;
>> 154 }
>> 155 else
>> 156 {
>> 157 theBuff.SetX(ii, GetX(high,i2));
>> 158 //*******************************************************************************
>> 159 //Change by E.Mendoza and D.Cano (CIEMAT):
>> 160 //y1 = GetY(high,i2);
>> 161 //y2 = GetY(low, GetX(high,i2)); //prob at ident theta
>> 162 y2 = GetY(high,i2);
>> 163 y1 = GetY(low, GetX(high,i2)); //prob at ident theta
>> 164 //*******************************************************************************
>> 165 theBuff.SetY(ii, theInt.Interpolate(theManager.GetScheme(high),
>> 166 x, x1, x2, y1, y2)); //energy interpol
>> 167 theBuff.SetScheme(ii, data[high].GetScheme(i2));
>> 168 if(std::abs(GetX(low,i1)-GetX(high,i2))<0.001) i1++;
>> 169 i2++;
>> 170 ii++;
>> 171 }
152 } 172 }
153 } 173 }
>> 174 //buff is full, now sample.
>> 175 return theBuff.Sample();
154 } 176 }
155 // buff is full, now sample. <<
156 return theBuff.Sample(); <<
157 } <<
158 177