| 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 // $Id: G4LivermoreComptonModel.cc 84806 2014-10-21 09:16:33Z gcosmo $
>> 27 // GEANT4 tag $Name: not supported by cvs2svn $
26 // 28 //
27 // 29 //
28 // Author: Alexander Bagulya 30 // Author: Alexander Bagulya
29 // 11 March 2012 31 // 11 March 2012
30 // on the base of G4LivermoreComptonMo 32 // on the base of G4LivermoreComptonModel
31 // 33 //
32 // History: 34 // History:
33 // -------- 35 // --------
34 // 18 Apr 2009 V Ivanchenko Cleanup initiali 36 // 18 Apr 2009 V Ivanchenko Cleanup initialisation and generation of secondaries:
35 // - apply internal high-ener << 37 // - apply internal high-energy limit only in constructor
36 // - do not apply low-energy 38 // - do not apply low-energy limit (default is 0)
37 // - remove GetMeanFreePath m 39 // - remove GetMeanFreePath method and table
38 // - added protection against << 40 // - added protection against numerical problem in energy sampling
39 // - use G4ElementSelector 41 // - use G4ElementSelector
40 // 26 Dec 2010 V Ivanchenko Load data tables 42 // 26 Dec 2010 V Ivanchenko Load data tables only once to avoid memory leak
41 43
42 #include "G4LivermoreComptonModel.hh" 44 #include "G4LivermoreComptonModel.hh"
43 <<
44 #include "G4AtomicShell.hh" <<
45 #include "G4AutoLock.hh" <<
46 #include "G4DopplerProfile.hh" <<
47 #include "G4Electron.hh" <<
48 #include "G4Exp.hh" <<
49 #include "G4Log.hh" <<
50 #include "G4LossTableManager.hh" <<
51 #include "G4ParticleChangeForGamma.hh" <<
52 #include "G4PhysicalConstants.hh" 45 #include "G4PhysicalConstants.hh"
53 #include "G4ShellData.hh" <<
54 #include "G4SystemOfUnits.hh" 46 #include "G4SystemOfUnits.hh"
>> 47 #include "G4Electron.hh"
>> 48 #include "G4ParticleChangeForGamma.hh"
>> 49 #include "G4LossTableManager.hh"
55 #include "G4VAtomDeexcitation.hh" 50 #include "G4VAtomDeexcitation.hh"
>> 51 #include "G4AtomicShell.hh"
>> 52 #include "G4Gamma.hh"
>> 53 #include "G4ShellData.hh"
>> 54 #include "G4DopplerProfile.hh"
>> 55 #include "G4Log.hh"
>> 56 #include "G4Exp.hh"
56 57
57 //....oooOO0OOooo........oooOO0OOooo........oo 58 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
58 59
59 namespace <<
60 { <<
61 G4Mutex LivermoreComptonModelMutex = G4MUTEX_I <<
62 } <<
63 using namespace std; 60 using namespace std;
64 61
65 G4PhysicsFreeVector* G4LivermoreComptonModel:: << 62 G4int G4LivermoreComptonModel::maxZ = 99;
66 G4ShellData* G4LivermoreComptonModel::shellDat << 63 G4LPhysicsFreeVector* G4LivermoreComptonModel::data[] = {0};
67 G4DopplerProfile* G4LivermoreComptonModel::pro << 64 G4ShellData* G4LivermoreComptonModel::shellData = 0;
68 G4String G4LivermoreComptonModel::gDataDirecto << 65 G4DopplerProfile* G4LivermoreComptonModel::profileData = 0;
69 66
70 static const G4double ln10 = G4Log(10.); 67 static const G4double ln10 = G4Log(10.);
71 68
72 G4LivermoreComptonModel::G4LivermoreComptonMod << 69 G4LivermoreComptonModel::G4LivermoreComptonModel(const G4ParticleDefinition*,
73 : G4VEmModel(nam) << 70 const G4String& nam)
74 { << 71 : G4VEmModel(nam),isInitialised(false)
75 fParticleChange = nullptr; << 72 {
76 verboseLevel = 1; << 73 lowestEnergy = 10 * eV;
>> 74
>> 75 verboseLevel=1 ;
77 // Verbosity scale: 76 // Verbosity scale:
78 // 0 = nothing << 77 // 0 = nothing
79 // 1 = warning for energy non-conservation << 78 // 1 = warning for energy non-conservation
80 // 2 = details of energy budget 79 // 2 = details of energy budget
81 // 3 = calculation of cross sections, file o 80 // 3 = calculation of cross sections, file openings, sampling of atoms
82 // 4 = entering in methods 81 // 4 = entering in methods
83 82
84 if (verboseLevel > 1) { << 83 if( verboseLevel>1 ) {
85 G4cout << "Livermore Compton model is cons 84 G4cout << "Livermore Compton model is constructed " << G4endl;
86 } 85 }
87 86
88 // Mark this model as "applicable" for atomi << 87 //Mark this model as "applicable" for atomic deexcitation
89 SetDeexcitationFlag(true); 88 SetDeexcitationFlag(true);
90 89
91 fParticleChange = nullptr; << 90 fParticleChange = 0;
92 fAtomDeexcitation = nullptr; << 91 fAtomDeexcitation = 0;
93 } 92 }
94 93
95 //....oooOO0OOooo........oooOO0OOooo........oo 94 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
96 95
97 G4LivermoreComptonModel::~G4LivermoreComptonMo 96 G4LivermoreComptonModel::~G4LivermoreComptonModel()
98 { << 97 {
99 if (IsMaster()) { << 98 if(IsMaster()) {
100 delete shellData; 99 delete shellData;
101 shellData = nullptr; << 100 shellData = 0;
102 delete profileData; 101 delete profileData;
103 profileData = nullptr; << 102 profileData = 0;
104 for (G4int i = 0; i <= maxZ; ++i) { <<
105 if (data[i]) { <<
106 delete data[i]; <<
107 data[i] = nullptr; <<
108 } <<
109 } <<
110 } 103 }
111 } 104 }
112 105
113 //....oooOO0OOooo........oooOO0OOooo........oo 106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
114 107
115 void G4LivermoreComptonModel::Initialise(const 108 void G4LivermoreComptonModel::Initialise(const G4ParticleDefinition* particle,
116 const << 109 const G4DataVector& cuts)
117 { 110 {
118 if (verboseLevel > 1) { 111 if (verboseLevel > 1) {
119 G4cout << "Calling G4LivermoreComptonModel 112 G4cout << "Calling G4LivermoreComptonModel::Initialise()" << G4endl;
120 } 113 }
121 114
122 // Initialise element selector 115 // Initialise element selector
123 if (IsMaster()) { << 116
124 // Initialise element selector << 117 if(IsMaster()) {
125 InitialiseElementSelectors(particle, cuts) <<
126 118
127 // Access to elements 119 // Access to elements
128 const G4ElementTable* elemTable = G4Elemen << 120
129 size_t numElems = (*elemTable).size(); << 121 char* path = getenv("G4LEDATA");
130 for (size_t ie = 0; ie < numElems; ++ie) { << 122
131 const G4Element* elem = (*elemTable)[ie] << 123 G4ProductionCutsTable* theCoupleTable =
132 const G4int Z = std::min(maxZ, elem->Get << 124 G4ProductionCutsTable::GetProductionCutsTable();
133 if (data[Z] == nullptr) { << 125 G4int numOfCouples = theCoupleTable->GetTableSize();
134 ReadData(Z); << 126
>> 127 for(G4int i=0; i<numOfCouples; ++i) {
>> 128 const G4Material* material =
>> 129 theCoupleTable->GetMaterialCutsCouple(i)->GetMaterial();
>> 130 const G4ElementVector* theElementVector = material->GetElementVector();
>> 131 G4int nelm = material->GetNumberOfElements();
>> 132
>> 133 for (G4int j=0; j<nelm; ++j) {
>> 134 G4int Z = G4lrint((*theElementVector)[j]->GetZ());
>> 135 if(Z < 1) { Z = 1; }
>> 136 else if(Z > maxZ){ Z = maxZ; }
>> 137
>> 138 if( (!data[Z]) ) { ReadData(Z, path); }
135 } 139 }
136 } 140 }
137 141
138 // For Doppler broadening 142 // For Doppler broadening
139 if (shellData == nullptr) { << 143 if(!shellData) {
140 shellData = new G4ShellData(); << 144 shellData = new G4ShellData();
141 shellData->SetOccupancyData(); 145 shellData->SetOccupancyData();
142 G4String file = "/doppler/shell-doppler" 146 G4String file = "/doppler/shell-doppler";
143 shellData->LoadData(file); 147 shellData->LoadData(file);
144 } 148 }
145 if (profileData == nullptr) { << 149 if(!profileData) { profileData = new G4DopplerProfile(); }
146 profileData = new G4DopplerProfile(); << 150
147 } << 151 InitialiseElementSelectors(particle, cuts);
148 } 152 }
149 153
150 if (verboseLevel > 2) { 154 if (verboseLevel > 2) {
151 G4cout << "Loaded cross section files" << 155 G4cout << "Loaded cross section files" << G4endl;
152 } 156 }
153 << 157
154 if (verboseLevel > 1) { << 158 if( verboseLevel>1 ) {
155 G4cout << "G4LivermoreComptonModel is init 159 G4cout << "G4LivermoreComptonModel is initialized " << G4endl
156 << "Energy range: " << LowEnergyLim << 160 << "Energy range: "
157 << " GeV" << G4endl; << 161 << LowEnergyLimit() / eV << " eV - "
158 } << 162 << HighEnergyLimit() / GeV << " GeV"
159 // << 163 << G4endl;
160 if (isInitialised) { <<
161 return; <<
162 } 164 }
>> 165 //
>> 166 if(isInitialised) { return; }
163 167
164 fParticleChange = GetParticleChangeForGamma( 168 fParticleChange = GetParticleChangeForGamma();
165 fAtomDeexcitation = G4LossTableManager::Inst << 169 fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
166 isInitialised = true; 170 isInitialised = true;
167 } 171 }
168 172
169 //....oooOO0OOooo........oooOO0OOooo........oo 173 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
170 174
171 void G4LivermoreComptonModel::InitialiseLocal( << 175 void G4LivermoreComptonModel::InitialiseLocal(const G4ParticleDefinition*,
>> 176 G4VEmModel* masterModel)
172 { 177 {
173 SetElementSelectors(masterModel->GetElementS 178 SetElementSelectors(masterModel->GetElementSelectors());
174 } 179 }
175 180
176 //....oooOO0OOooo........oooOO0OOooo........oo 181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
177 182
178 const G4String& G4LivermoreComptonModel::FindD << 183 void G4LivermoreComptonModel::ReadData(size_t Z, const char* path)
179 { <<
180 // no check in this method - environment var <<
181 if (gDataDirectory.empty()) { <<
182 auto param = G4EmParameters::Instance(); <<
183 std::ostringstream ost; <<
184 if (param->LivermoreDataDir() == "livermor <<
185 ost << param->GetDirLEDATA() << "/liverm <<
186 } <<
187 else { <<
188 ost << param->GetDirLEDATA() << "/epics2 <<
189 } <<
190 gDataDirectory = ost.str(); <<
191 } <<
192 return gDataDirectory; <<
193 } <<
194 <<
195 //....oooOO0OOooo........oooOO0OOooo........oo <<
196 <<
197 void G4LivermoreComptonModel::ReadData(const G <<
198 { 184 {
199 if (verboseLevel > 1) { << 185 if (verboseLevel > 1)
200 G4cout << "G4LivermoreComptonModel::ReadDa << 186 {
201 } << 187 G4cout << "G4LivermoreComptonModel::ReadData()"
202 const G4int Z = std::min(ZZ, maxZ); << 188 << G4endl;
203 << 189 }
204 if (data[Z] != nullptr) { << 190 if(data[Z]) { return; }
205 return; << 191 const char* datadir = path;
206 } << 192 if(!datadir)
207 << 193 {
208 data[Z] = new G4PhysicsFreeVector(); << 194 datadir = getenv("G4LEDATA");
209 << 195 if(!datadir)
>> 196 {
>> 197 G4Exception("G4LivermoreComptonModel::ReadData()",
>> 198 "em0006",FatalException,
>> 199 "Environment variable G4LEDATA not defined");
>> 200 return;
>> 201 }
>> 202 }
>> 203
>> 204 data[Z] = new G4LPhysicsFreeVector();
>> 205
>> 206 // Activation of spline interpolation
>> 207 data[Z]->SetSpline(false);
>> 208
210 std::ostringstream ost; 209 std::ostringstream ost;
211 ost << FindDirectoryPath() << "ce-cs-" << Z << 210 ost << datadir << "/livermore/comp/ce-cs-" << Z <<".dat";
212 <<
213 std::ifstream fin(ost.str().c_str()); 211 std::ifstream fin(ost.str().c_str());
214 << 212
215 if (!fin.is_open()) { << 213 if( !fin.is_open())
216 G4ExceptionDescription ed; << 214 {
217 ed << "G4LivermoreComptonModel data file < << 215 G4ExceptionDescription ed;
218 << G4endl; << 216 ed << "G4LivermoreComptonModel data file <" << ost.str().c_str()
219 G4Exception("G4LivermoreComptonModel::Read << 217 << "> is not opened!" << G4endl;
220 "G4LEDATA version should be G4 << 218 G4Exception("G4LivermoreComptonModel::ReadData()",
221 return; << 219 "em0003",FatalException,
222 } << 220 ed,"G4LEDATA version should be G4EMLOW6.34 or later");
223 else { << 221 return;
224 if (verboseLevel > 3) { << 222 } else {
225 G4cout << "File " << ost.str() << " is o << 223 if(verboseLevel > 3) {
226 } << 224 G4cout << "File " << ost.str()
227 data[Z]->Retrieve(fin, true); << 225 << " is opened by G4LivermoreComptonModel" << G4endl;
228 data[Z]->ScaleVector(MeV, MeV * barn); << 226 }
229 } << 227 data[Z]->Retrieve(fin, true);
>> 228 data[Z]->ScaleVector(MeV, MeV*barn);
>> 229 }
230 fin.close(); 230 fin.close();
231 } 231 }
232 232
233 //....oooOO0OOooo........oooOO0OOooo........oo 233 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
234 234
235 G4double G4LivermoreComptonModel::ComputeCross << 235 G4double
236 << 236 G4LivermoreComptonModel::ComputeCrossSectionPerAtom(const G4ParticleDefinition*,
237 << 237 G4double GammaEnergy,
>> 238 G4double Z, G4double,
>> 239 G4double, G4double)
238 { 240 {
239 if (verboseLevel > 3) { 241 if (verboseLevel > 3) {
240 G4cout << "G4LivermoreComptonModel::Comput << 242 G4cout << "G4LivermoreComptonModel::ComputeCrossSectionPerAtom()"
>> 243 << G4endl;
241 } 244 }
242 G4double cs = 0.0; << 245 G4double cs = 0.0;
243 246
244 if (GammaEnergy < LowEnergyLimit()) { << 247 if (GammaEnergy < lowestEnergy) { return 0.0; }
245 return 0.0; <<
246 } <<
247 248
248 G4int intZ = G4lrint(Z); 249 G4int intZ = G4lrint(Z);
249 if (intZ < 1 || intZ > maxZ) { << 250 if(intZ < 1 || intZ > maxZ) { return cs; }
250 return cs; <<
251 } <<
252 251
253 G4PhysicsFreeVector* pv = data[intZ]; << 252 G4LPhysicsFreeVector* pv = data[intZ];
254 253
255 // if element was not initialised 254 // if element was not initialised
256 // do initialisation safely for MT mode 255 // do initialisation safely for MT mode
257 if (pv == nullptr) { << 256 if(!pv)
258 InitialiseForElement(nullptr, intZ); << 257 {
259 pv = data[intZ]; << 258 InitialiseForElement(0, intZ);
260 if (pv == nullptr) { << 259 pv = data[intZ];
261 return cs; << 260 if(!pv) { return cs; }
262 } 261 }
263 } <<
264 262
265 auto n = G4int(pv->GetVectorLength() - 1); << 263 G4int n = pv->GetVectorLength() - 1;
266 G4double e1 = pv->Energy(0); 264 G4double e1 = pv->Energy(0);
267 G4double e2 = pv->Energy(n); 265 G4double e2 = pv->Energy(n);
268 266
269 if (GammaEnergy <= e1) { << 267 if(GammaEnergy <= e1) { cs = GammaEnergy/(e1*e1)*pv->Value(e1); }
270 cs = GammaEnergy / (e1 * e1) * pv->Value(e << 268 else if(GammaEnergy <= e2) { cs = pv->Value(GammaEnergy)/GammaEnergy; }
271 } << 269 else if(GammaEnergy > e2) { cs = pv->Value(e2)/GammaEnergy; }
272 else if (GammaEnergy <= e2) { << 270
273 cs = pv->Value(GammaEnergy) / GammaEnergy; <<
274 } <<
275 else if (GammaEnergy > e2) { <<
276 cs = pv->Value(e2) / GammaEnergy; <<
277 } <<
278 <<
279 return cs; 271 return cs;
280 } 272 }
281 273
282 //....oooOO0OOooo........oooOO0OOooo........oo 274 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
283 275
284 void G4LivermoreComptonModel::SampleSecondarie << 276
285 << 277 void G4LivermoreComptonModel::SampleSecondaries(
286 << 278 std::vector<G4DynamicParticle*>* fvect,
287 << 279 const G4MaterialCutsCouple* couple,
>> 280 const G4DynamicParticle* aDynamicGamma,
>> 281 G4double, G4double)
288 { 282 {
289 // The scattered gamma energy is sampled acc << 283
290 // formula then accepted or rejected dependi << 284 // The scattered gamma energy is sampled according to Klein - Nishina
>> 285 // formula then accepted or rejected depending on the Scattering Function
291 // multiplied by factor from Klein - Nishina 286 // multiplied by factor from Klein - Nishina formula.
292 // Expression of the angular distribution as 287 // Expression of the angular distribution as Klein Nishina
293 // angular and energy distribution and Scatt 288 // angular and energy distribution and Scattering fuctions is taken from
294 // D. E. Cullen "A simple model of photon tr 289 // D. E. Cullen "A simple model of photon transport" Nucl. Instr. Meth.
295 // Phys. Res. B 101 (1995). Method of sampli 290 // Phys. Res. B 101 (1995). Method of sampling with form factors is different
296 // data are interpolated while in the articl 291 // data are interpolated while in the article they are fitted.
297 // Reference to the article is from J. Stepa 292 // Reference to the article is from J. Stepanek New Photon, Positron
298 // and Electron Interaction Data for GEANT i 293 // and Electron Interaction Data for GEANT in Energy Range from 1 eV to 10
299 // TeV (draft). 294 // TeV (draft).
300 // The random number techniques of Butcher & 295 // The random number techniques of Butcher & Messel are used
301 // (Nucl Phys 20(1960),15). 296 // (Nucl Phys 20(1960),15).
302 297
303 G4double photonEnergy0 = aDynamicGamma->GetK 298 G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
304 299
305 if (verboseLevel > 3) { 300 if (verboseLevel > 3) {
306 G4cout << "G4LivermoreComptonModel::Sample << 301 G4cout << "G4LivermoreComptonModel::SampleSecondaries() E(MeV)= "
307 << " in " << couple->GetMaterial()- << 302 << photonEnergy0/MeV << " in " << couple->GetMaterial()->GetName()
>> 303 << G4endl;
308 } 304 }
309 305
310 // do nothing below the threshold << 306 // low-energy gamma is absorpted by this process
311 // should never get here because the XS is z << 307 if (photonEnergy0 <= lowestEnergy)
312 if (photonEnergy0 < LowEnergyLimit()) return << 308 {
>> 309 fParticleChange->ProposeTrackStatus(fStopAndKill);
>> 310 fParticleChange->SetProposedKineticEnergy(0.);
>> 311 fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
>> 312 return ;
>> 313 }
313 314
314 G4double e0m = photonEnergy0 / electron_mass << 315 G4double e0m = photonEnergy0 / electron_mass_c2 ;
315 G4ParticleMomentum photonDirection0 = aDynam 316 G4ParticleMomentum photonDirection0 = aDynamicGamma->GetMomentumDirection();
316 317
317 // Select randomly one element in the curren 318 // Select randomly one element in the current material
318 const G4ParticleDefinition* particle = aDyna << 319 const G4ParticleDefinition* particle = aDynamicGamma->GetDefinition();
319 const G4Element* elm = SelectRandomAtom(coup << 320 const G4Element* elm = SelectRandomAtom(couple,particle,photonEnergy0);
320 321
321 G4int Z = elm->GetZasInt(); << 322 G4int Z = G4lrint(elm->GetZ());
322 323
323 G4double epsilon0Local = 1. / (1. + 2. * e0m 324 G4double epsilon0Local = 1. / (1. + 2. * e0m);
324 G4double epsilon0Sq = epsilon0Local * epsilo 325 G4double epsilon0Sq = epsilon0Local * epsilon0Local;
325 G4double alpha1 = -G4Log(epsilon0Local); 326 G4double alpha1 = -G4Log(epsilon0Local);
326 G4double alpha2 = 0.5 * (1. - epsilon0Sq); 327 G4double alpha2 = 0.5 * (1. - epsilon0Sq);
327 328
328 G4double wlPhoton = h_Planck * c_light / pho << 329 G4double wlPhoton = h_Planck*c_light/photonEnergy0;
329 330
330 // Sample the energy of the scattered photon 331 // Sample the energy of the scattered photon
331 G4double epsilon; 332 G4double epsilon;
332 G4double epsilonSq; 333 G4double epsilonSq;
333 G4double oneCosT; 334 G4double oneCosT;
334 G4double sinT2; 335 G4double sinT2;
335 G4double gReject; 336 G4double gReject;
336 337
337 if (verboseLevel > 3) { 338 if (verboseLevel > 3) {
338 G4cout << "Started loop to sample gamma en 339 G4cout << "Started loop to sample gamma energy" << G4endl;
339 } 340 }
340 << 341
341 do { 342 do {
342 if (alpha1 / (alpha1 + alpha2) > G4Uniform << 343 if ( alpha1/(alpha1+alpha2) > G4UniformRand())
343 epsilon = G4Exp(-alpha1 * G4UniformRand( << 344 {
344 epsilonSq = epsilon * epsilon; << 345 epsilon = G4Exp(-alpha1 * G4UniformRand());
345 } << 346 epsilonSq = epsilon * epsilon;
346 else { << 347 }
347 epsilonSq = epsilon0Sq + (1. - epsilon0S << 348 else
348 epsilon = std::sqrt(epsilonSq); << 349 {
349 } << 350 epsilonSq = epsilon0Sq + (1. - epsilon0Sq) * G4UniformRand();
>> 351 epsilon = std::sqrt(epsilonSq);
>> 352 }
350 353
351 oneCosT = (1. - epsilon) / (epsilon * e0m) << 354 oneCosT = (1. - epsilon) / ( epsilon * e0m);
352 sinT2 = oneCosT * (2. - oneCosT); << 355 sinT2 = oneCosT * (2. - oneCosT);
353 G4double x = std::sqrt(oneCosT / 2.) * cm << 356 G4double x = std::sqrt(oneCosT/2.) * cm / wlPhoton;
354 G4double scatteringFunction = ComputeScatt << 357 G4double scatteringFunction = ComputeScatteringFunction(x, Z);
355 gReject = (1. - epsilon * sinT2 / (1. + ep << 358 gReject = (1. - epsilon * sinT2 / (1. + epsilonSq)) * scatteringFunction;
356 359
357 } while (gReject < G4UniformRand() * Z); << 360 } while(gReject < G4UniformRand()*Z);
358 361
359 G4double cosTheta = 1. - oneCosT; 362 G4double cosTheta = 1. - oneCosT;
360 G4double sinTheta = std::sqrt(sinT2); << 363 G4double sinTheta = std::sqrt (sinT2);
361 G4double phi = twopi * G4UniformRand(); << 364 G4double phi = twopi * G4UniformRand() ;
362 G4double dirx = sinTheta * std::cos(phi); 365 G4double dirx = sinTheta * std::cos(phi);
363 G4double diry = sinTheta * std::sin(phi); 366 G4double diry = sinTheta * std::sin(phi);
364 G4double dirz = cosTheta; << 367 G4double dirz = cosTheta ;
365 368
366 // Doppler broadening - Method based on: 369 // Doppler broadening - Method based on:
367 // Y. Namito, S. Ban and H. Hirayama, << 370 // Y. Namito, S. Ban and H. Hirayama,
368 // "Implementation of the Doppler Broadening << 371 // "Implementation of the Doppler Broadening of a Compton-Scattered Photon
369 // into the EGS4 Code", NIM A 349, pp. 489-4 372 // into the EGS4 Code", NIM A 349, pp. 489-494, 1994
370 << 373
371 // Maximum number of sampling iterations 374 // Maximum number of sampling iterations
372 static G4int maxDopplerIterations = 1000; 375 static G4int maxDopplerIterations = 1000;
373 G4double bindingE = 0.; 376 G4double bindingE = 0.;
374 G4double photonEoriginal = epsilon * photonE 377 G4double photonEoriginal = epsilon * photonEnergy0;
375 G4double photonE = -1.; 378 G4double photonE = -1.;
376 G4int iteration = 0; 379 G4int iteration = 0;
377 G4double eMax = photonEnergy0; 380 G4double eMax = photonEnergy0;
378 381
379 G4int shellIdx = 0; 382 G4int shellIdx = 0;
380 383
381 if (verboseLevel > 3) { 384 if (verboseLevel > 3) {
382 G4cout << "Started loop to sample broading 385 G4cout << "Started loop to sample broading" << G4endl;
383 } 386 }
384 387
385 do { << 388 do
386 ++iteration; << 389 {
387 // Select shell based on shell occupancy << 390 ++iteration;
388 shellIdx = shellData->SelectRandomShell(Z) << 391 // Select shell based on shell occupancy
389 bindingE = shellData->BindingEnergy(Z, she << 392 shellIdx = shellData->SelectRandomShell(Z);
390 << 393 bindingE = shellData->BindingEnergy(Z,shellIdx);
391 if (verboseLevel > 3) { << 394
392 G4cout << "Shell ID= " << shellIdx << " << 395 if (verboseLevel > 3) {
393 } << 396 G4cout << "Shell ID= " << shellIdx
394 << 397 << " Ebind(keV)= " << bindingE/keV << G4endl;
395 eMax = photonEnergy0 - bindingE; <<
396 <<
397 // Randomly sample bound electron momentum <<
398 // (memento: the data set is in Atomic Uni <<
399 G4double pSample = profileData->RandomSele <<
400 if (verboseLevel > 3) { <<
401 G4cout << "pSample= " << pSample << G4en <<
402 } <<
403 // Rescale from atomic units <<
404 G4double pDoppler = pSample * fine_structu <<
405 G4double pDoppler2 = pDoppler * pDoppler; <<
406 G4double var2 = 1. + oneCosT * e0m; <<
407 G4double var3 = var2 * var2 - pDoppler2; <<
408 G4double var4 = var2 - pDoppler2 * cosThet <<
409 G4double var = var4 * var4 - var3 + pDoppl <<
410 if (var > 0.) { <<
411 G4double varSqrt = std::sqrt(var); <<
412 G4double scale = photonEnergy0 / var3; <<
413 // Random select either root <<
414 if (G4UniformRand() < 0.5) { <<
415 photonE = (var4 - varSqrt) * scale; <<
416 } 398 }
417 else { << 399
418 photonE = (var4 + varSqrt) * scale; << 400 eMax = photonEnergy0 - bindingE;
>> 401
>> 402 // Randomly sample bound electron momentum
>> 403 // (memento: the data set is in Atomic Units)
>> 404 G4double pSample = profileData->RandomSelectMomentum(Z,shellIdx);
>> 405 if (verboseLevel > 3) {
>> 406 G4cout << "pSample= " << pSample << G4endl;
419 } 407 }
420 } << 408 // Rescale from atomic units
421 else { << 409 G4double pDoppler = pSample * fine_structure_const;
422 photonE = -1.; << 410 G4double pDoppler2 = pDoppler * pDoppler;
423 } << 411 G4double var2 = 1. + oneCosT * e0m;
424 } while (iteration <= maxDopplerIterations & << 412 G4double var3 = var2*var2 - pDoppler2;
425 << 413 G4double var4 = var2 - pDoppler2 * cosTheta;
>> 414 G4double var = var4*var4 - var3 + pDoppler2 * var3;
>> 415 if (var > 0.)
>> 416 {
>> 417 G4double varSqrt = std::sqrt(var);
>> 418 G4double scale = photonEnergy0 / var3;
>> 419 // Random select either root
>> 420 if (G4UniformRand() < 0.5) { photonE = (var4 - varSqrt) * scale; }
>> 421 else { photonE = (var4 + varSqrt) * scale; }
>> 422 }
>> 423 else
>> 424 {
>> 425 photonE = -1.;
>> 426 }
>> 427 } while (iteration <= maxDopplerIterations && photonE > eMax);
>> 428 // (photonE < 0. || photonE > eMax || photonE < eMax*G4UniformRand()) );
>> 429
>> 430
426 // End of recalculation of photon energy wit 431 // End of recalculation of photon energy with Doppler broadening
427 // Revert to original if maximum number of i << 432 // Revert to original if maximum number of iterations threshold
428 // has been reached 433 // has been reached
429 if (iteration >= maxDopplerIterations) { << 434
430 photonE = photonEoriginal; << 435 if (iteration >= maxDopplerIterations)
431 bindingE = 0.; << 436 {
432 } << 437 photonE = photonEoriginal;
>> 438 bindingE = 0.;
>> 439 }
433 440
434 // Update G4VParticleChange for the scattere 441 // Update G4VParticleChange for the scattered photon
435 G4ThreeVector photonDirection1(dirx, diry, d << 442
>> 443 G4ThreeVector photonDirection1(dirx,diry,dirz);
436 photonDirection1.rotateUz(photonDirection0); 444 photonDirection1.rotateUz(photonDirection0);
437 fParticleChange->ProposeMomentumDirection(ph << 445 fParticleChange->ProposeMomentumDirection(photonDirection1) ;
438 446
439 G4double photonEnergy1 = photonE; 447 G4double photonEnergy1 = photonE;
440 448
441 if (photonEnergy1 > 0.) { 449 if (photonEnergy1 > 0.) {
442 fParticleChange->SetProposedKineticEnergy( << 450 fParticleChange->SetProposedKineticEnergy(photonEnergy1) ;
443 } << 451
444 else { << 452 } else {
445 // photon absorbed 453 // photon absorbed
446 photonEnergy1 = 0.; 454 photonEnergy1 = 0.;
447 fParticleChange->SetProposedKineticEnergy( << 455 fParticleChange->SetProposedKineticEnergy(0.) ;
448 fParticleChange->ProposeTrackStatus(fStopA << 456 fParticleChange->ProposeTrackStatus(fStopAndKill);
449 fParticleChange->ProposeLocalEnergyDeposit 457 fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
450 return; 458 return;
451 } 459 }
452 460
453 // Kinematics of the scattered electron 461 // Kinematics of the scattered electron
454 G4double eKineticEnergy = photonEnergy0 - ph 462 G4double eKineticEnergy = photonEnergy0 - photonEnergy1 - bindingE;
455 463
456 // protection against negative final energy: 464 // protection against negative final energy: no e- is created
457 if (eKineticEnergy < 0.0) { << 465 if(eKineticEnergy < 0.0) {
458 fParticleChange->ProposeLocalEnergyDeposit 466 fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0 - photonEnergy1);
459 return; 467 return;
460 } 468 }
461 469
462 G4double eTotalEnergy = eKineticEnergy + ele 470 G4double eTotalEnergy = eKineticEnergy + electron_mass_c2;
463 471
464 G4double electronE = photonEnergy0 * (1. - e << 472 G4double electronE = photonEnergy0 * (1. - epsilon) + electron_mass_c2;
465 G4double electronP2 = electronE * electronE << 473 G4double electronP2 =
>> 474 electronE*electronE - electron_mass_c2*electron_mass_c2;
466 G4double sinThetaE = -1.; 475 G4double sinThetaE = -1.;
467 G4double cosThetaE = 0.; 476 G4double cosThetaE = 0.;
468 if (electronP2 > 0.) { << 477 if (electronP2 > 0.)
469 cosThetaE = (eTotalEnergy + photonEnergy1) << 478 {
470 sinThetaE = -1. * sqrt(1. - cosThetaE * co << 479 cosThetaE = (eTotalEnergy + photonEnergy1 )*
471 } << 480 (1. - epsilon) / std::sqrt(electronP2);
472 << 481 sinThetaE = -1. * sqrt(1. - cosThetaE * cosThetaE);
>> 482 }
>> 483
473 G4double eDirX = sinThetaE * std::cos(phi); 484 G4double eDirX = sinThetaE * std::cos(phi);
474 G4double eDirY = sinThetaE * std::sin(phi); 485 G4double eDirY = sinThetaE * std::sin(phi);
475 G4double eDirZ = cosThetaE; 486 G4double eDirZ = cosThetaE;
476 487
477 G4ThreeVector eDirection(eDirX, eDirY, eDirZ << 488 G4ThreeVector eDirection(eDirX,eDirY,eDirZ);
478 eDirection.rotateUz(photonDirection0); 489 eDirection.rotateUz(photonDirection0);
479 auto dp = new G4DynamicParticle(G4Electron:: << 490 G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),
>> 491 eDirection,eKineticEnergy) ;
480 fvect->push_back(dp); 492 fvect->push_back(dp);
481 493
482 // sample deexcitation 494 // sample deexcitation
>> 495 //
>> 496
483 if (verboseLevel > 3) { 497 if (verboseLevel > 3) {
484 G4cout << "Started atomic de-excitation " 498 G4cout << "Started atomic de-excitation " << fAtomDeexcitation << G4endl;
485 } 499 }
486 500
487 if (nullptr != fAtomDeexcitation && iteratio << 501 if(fAtomDeexcitation && iteration < maxDopplerIterations) {
488 G4int index = couple->GetIndex(); 502 G4int index = couple->GetIndex();
489 if (fAtomDeexcitation->CheckDeexcitationAc << 503 if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
490 size_t nbefore = fvect->size(); 504 size_t nbefore = fvect->size();
491 auto as = G4AtomicShellEnumerator(shellI << 505 G4AtomicShellEnumerator as = G4AtomicShellEnumerator(shellIdx);
492 const G4AtomicShell* shell = fAtomDeexci 506 const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
493 fAtomDeexcitation->GenerateParticles(fve 507 fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
494 size_t nafter = fvect->size(); 508 size_t nafter = fvect->size();
495 if (nafter > nbefore) { << 509 if(nafter > nbefore) {
496 for (size_t i = nbefore; i < nafter; + << 510 for (size_t i=nbefore; i<nafter; ++i) {
497 // Check if there is enough residual << 511 //Check if there is enough residual energy
498 if (bindingE >= ((*fvect)[i])->GetKi << 512 if (bindingE >= ((*fvect)[i])->GetKineticEnergy())
499 // Ok, this is a valid secondary: << 513 {
500 bindingE -= ((*fvect)[i])->GetKine << 514 //Ok, this is a valid secondary: keep it
501 } << 515 bindingE -= ((*fvect)[i])->GetKineticEnergy();
502 else { << 516 }
503 // Invalid secondary: not enough e << 517 else
504 // Keep its energy in the local de << 518 {
505 delete (*fvect)[i]; << 519 //Invalid secondary: not enough energy to create it!
506 (*fvect)[i] = nullptr; << 520 //Keep its energy in the local deposit
507 } << 521 delete (*fvect)[i];
508 } << 522 (*fvect)[i]=0;
>> 523 }
>> 524 }
509 } 525 }
510 } 526 }
511 } 527 }
512 bindingE = std::max(bindingE, 0.0); << 528
>> 529 //This should never happen
>> 530 if(bindingE < 0.0)
>> 531 G4Exception("G4LowEPComptonModel::SampleSecondaries()",
>> 532 "em2051",FatalException,"Negative local energy deposit");
>> 533
513 fParticleChange->ProposeLocalEnergyDeposit(b 534 fParticleChange->ProposeLocalEnergyDeposit(bindingE);
>> 535
514 } 536 }
515 537
516 //....oooOO0OOooo........oooOO0OOooo........oo 538 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
517 539
518 G4double G4LivermoreComptonModel::ComputeScatt << 540 G4double
>> 541 G4LivermoreComptonModel::ComputeScatteringFunction(G4double x, G4int Z)
519 { 542 {
520 G4double value = Z; 543 G4double value = Z;
521 if (x <= ScatFuncFitParam[Z][3]) { << 544 if (x <= ScatFuncFitParam[Z][2]) {
522 G4double lgq = G4Log(x) / ln10; <<
523 545
524 if (lgq < ScatFuncFitParam[Z][1]) { << 546 G4double lgq = G4Log(x)/ln10;
525 value = ScatFuncFitParam[Z][4] + lgq * S << 547
526 } << 548 if (lgq < ScatFuncFitParam[Z][1]) {
527 else if (lgq >= ScatFuncFitParam[Z][1] && << 549 value = ScatFuncFitParam[Z][3] + lgq*ScatFuncFitParam[Z][4];
528 value = ScatFuncFitParam[Z][6] << 550 } else {
529 + lgq << 551 value = ScatFuncFitParam[Z][5] + lgq*ScatFuncFitParam[Z][6] +
530 * (ScatFuncFitParam[Z][7] << 552 lgq*lgq*ScatFuncFitParam[Z][7] + lgq*lgq*lgq*ScatFuncFitParam[Z][8];
531 + lgq <<
532 * (ScatFuncFitParam[Z <<
533 + lgq * (ScatFuncF <<
534 } 553 }
535 else { << 554 value = G4Exp(value*ln10);
536 value = ScatFuncFitParam[Z][11] <<
537 + lgq <<
538 * (ScatFuncFitParam[Z][12] <<
539 + lgq <<
540 * (ScatFuncFitParam[Z <<
541 + lgq * (ScatFuncF <<
542 } <<
543 value = G4Exp(value * ln10); <<
544 } 555 }
545 // G4cout << " value= " << value << G4end << 556 //G4cout << " value= " << value << G4endl;
546 return value; 557 return value;
547 } 558 }
548 559
549 //....oooOO0OOooo........oooOO0OOooo........oo 560 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 561 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 562
>> 563 #include "G4AutoLock.hh"
>> 564 namespace { G4Mutex LivermoreComptonModelMutex = G4MUTEX_INITIALIZER; }
550 565
551 void G4LivermoreComptonModel::InitialiseForEle << 566 void
>> 567 G4LivermoreComptonModel::InitialiseForElement(const G4ParticleDefinition*,
>> 568 G4int Z)
552 { 569 {
553 G4AutoLock l(&LivermoreComptonModelMutex); 570 G4AutoLock l(&LivermoreComptonModelMutex);
554 if (data[Z] == nullptr) { << 571 // G4cout << "G4LivermoreComptonModel::InitialiseForElement Z= "
555 ReadData(Z); << 572 // << Z << G4endl;
556 } << 573 if(!data[Z]) { ReadData(Z); }
557 l.unlock(); 574 l.unlock();
558 } 575 }
559 576
560 //....oooOO0OOooo........oooOO0OOooo........oo 577 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
561 578
562 // Fitting data to compute scattering function << 579 //Fitting data to compute scattering function
563 const G4double G4LivermoreComptonModel::ScatFu << 580 const G4double G4LivermoreComptonModel::ScatFuncFitParam[101][9] = {
564 {0, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., << 581 { 0, 0., 0., 0., 0., 0., 0., 0., 0.},
565 {1, 6.000000000e+00, 7.087999300e+00, 1.4996 << 582 { 1, 6.673, 1.49968E+08, -14.352, 1.999, -143.374, 50.787, -5.951, 0.2304418},
566 -3.925518125e+02, 2.434944521e+02, -5.78439 << 583 { 2, 6.500, 2.50035E+08, -14.215, 1.970, -53.649, 13.892, -0.948, 0.006996759},
567 -1.649463594e+03, 8.121933215e+02, -1.49831 << 584 { 3, 6.551, 3.99945E+08, -13.555, 1.993, -62.090, 21.462, -2.453, 0.093416},
568 {2, 6.000000000e+00, 7.199000403e+00, 2.5003 << 585 { 4, 6.500, 5.00035E+08, -13.746, 1.998, -127.906, 46.491, -5.614, 0.2262103},
569 3.574019365e+02, -1.978574937e+02, 3.971327 << 586 { 5, 6.500, 5.99791E+08, -13.800, 1.998, -131.153, 47.132, -5.619, 0.2233819},
570 -4.009960832e+02, 1.575831469e+02, -2.17476 << 587 { 6, 6.708, 6.99842E+08, -13.885, 1.999, -128.143, 45.379, -5.325, 0.2083009},
571 {3, 6.000000000e+00, 7.301000136e+00, 3.9994 << 588 { 7, 6.685, 7.99834E+08, -13.885, 2.000, -131.048, 46.314, -5.421, 0.2114925},
572 7.051635443e+02, -4.223841786e+02, 9.318729 << 589 { 8, 6.669, 7.99834E+08, -13.962, 2.001, -128.225, 44.818, -5.183, 0.1997155},
573 1.524679907e+03, -7.851479582e+02, 1.509941 << 590 { 9, 6.711, 7.99834E+08, -13.999, 2.000, -122.112, 42.103, -4.796, 0.1819099},
574 {4, 6.000000000e+00, 7.349500202e+00, 5.0003 << 591 { 10, 6.702, 7.99834E+08, -14.044, 1.999, -110.143, 37.225, -4.143, 0.1532094},
575 -1.832909604e+02, 1.193997722e+02, -3.03432 << 592 { 11, 6.425, 1.00000E+09, -13.423, 1.993, -41.137, 12.313, -1.152, 0.03384553},
576 1.397476657e+03, -7.026416933e+02, 1.320720 << 593 { 12, 6.542, 1.00000E+09, -13.389, 1.997, -53.549, 17.420, -1.840, 0.06431849},
577 {5, 6.000000000e+00, 7.388999972e+00, 5.9979 << 594 { 13, 6.570, 1.49968E+09, -13.401, 1.997, -66.243, 22.297, -2.460, 0.09045854},
578 -2.334197545e+02, 1.467013466e+02, -3.57485 << 595 { 14, 6.364, 1.49968E+09, -13.452, 1.999, -78.271, 26.757, -3.008, 0.1128195},
579 6.784713308e+02, -3.419562074e+02, 6.433945 << 596 { 15, 6.500, 1.49968E+09, -13.488, 1.998, -85.069, 29.164, -3.291, 0.1239113},
580 {6, 6.000000000e+00, 7.422500001e+00, 6.9984 << 597 { 16, 6.500, 1.49968E+09, -13.532, 1.998, -93.640, 32.274, -3.665, 0.1388633},
581 -2.460254935e+02, 1.516613633e+02, -3.62202 << 598 { 17, 6.500, 1.49968E+09, -13.584, 2.000, -98.534, 33.958, -3.857, 0.1461557},
582 -1.610185428e+02, 7.010907070e+01, -1.14237 << 599 { 18, 6.500, 1.49968E+09, -13.618, 1.999, -100.077, 34.379, -3.891, 0.1468902},
583 {7, 6.000000000e+00, 7.451499931e+00, 7.9983 << 600 { 19, 6.500, 1.99986E+09, -13.185, 1.992, -53.819, 17.528, -1.851, 0.0648722},
584 -3.054540719e+02, 1.877740247e+02, -4.44027 << 601 { 20, 6.490, 1.99986E+09, -13.123, 1.993, -52.221, 17.169, -1.832, 0.06502094},
585 -2.263864349e+02, 1.017885461e+02, -1.71698 << 602 { 21, 6.498, 1.99986E+09, -13.157, 1.994, -55.365, 18.276, -1.961, 0.07002778},
586 {8, 6.000000000e+00, 7.451499931e+00, 7.9983 << 603 { 22, 6.495, 1.99986E+09, -13.183, 1.994, -57.412, 18.957, -2.036, 0.07278856},
587 -3.877174895e+02, 2.345831969e+02, -5.43182 << 604 { 23, 6.487, 1.99986E+09, -13.216, 1.995, -58.478, 19.270, -2.065, 0.07362722},
588 -7.949384302e+02, 3.757293602e+02, -6.66174 << 605 { 24, 6.500, 1.99986E+09, -13.330, 1.997, -62.192, 20.358, -2.167, 0.07666583},
589 {9, 6.000000000e+00, 7.451499931e+00, 7.9983 << 606 { 25, 6.488, 1.99986E+09, -13.277, 1.997, -58.007, 18.924, -2.003, 0.0704305},
590 -2.939854827e+02, 1.784214589e+02, -4.16847 << 607 { 26, 6.500, 5.00035E+09, -13.292, 1.997, -61.176, 20.067, -2.141, 0.0760269},
591 -1.169326170e+03, 5.545642014e+02, -9.86302 << 608 { 27, 6.500, 5.00035E+09, -13.321, 1.998, -61.909, 20.271, -2.159, 0.07653559},
592 {10, 6.000000000e+00, 7.451499931e+00, 7.998 << 609 { 28, 6.500, 5.00035E+09, -13.340, 1.998, -62.402, 20.391, -2.167, 0.07664243},
593 -2.615701853e+02, 1.582596311e+02, -3.69811 << 610 { 29, 6.500, 5.00035E+09, -13.439, 1.998, -67.305, 21.954, -2.331, 0.0823267},
594 -1.275287356e+03, 6.022076554e+02, -1.06641 << 611 { 30, 6.500, 5.00035E+09, -13.383, 1.999, -62.064, 20.136, -2.122, 0.07437589},
595 {11, 6.000000000e+00, 7.500000000e+00, 1.000 << 612 { 31, 6.500, 5.00035E+09, -13.349, 1.997, -61.068, 19.808, -2.086, 0.07307488},
596 1.112662501e+03, -6.807056448e+02, 1.545837 << 613 { 32, 6.500, 5.00035E+09, -13.373, 1.999, -63.126, 20.553, -2.175, 0.07660222},
597 -1.007702307e+03, 4.699937040e+02, -8.22035 << 614 { 33, 6.500, 5.00035E+09, -13.395, 1.999, -65.674, 21.445, -2.278, 0.08054694},
598 {12, 6.000000000e+00, 7.500000000e+00, 1.000 << 615 { 34, 6.500, 5.00035E+09, -13.417, 1.999, -69.457, 22.811, -2.442, 0.08709536},
599 9.895649717e+02, -5.983228286e+02, 1.340681 << 616 { 35, 6.500, 5.00035E+09, -13.442, 2.000, -72.283, 23.808, -2.558, 0.09156808},
600 -5.790532602e+02, 2.626052403e+02, -4.46354 << 617 { 36, 6.500, 5.00035E+09, -13.451, 1.998, -74.696, 24.641, -2.653, 0.09516597},
601 {13, 6.000000000e+00, 7.587999300e+00, 1.499 << 618 { 37, 6.500, 5.00035E+09, -13.082, 1.991, -46.235, 14.519, -1.458, 0.04837659},
602 7.335256091e+02, -4.405291562e+02, 9.770954 << 619 { 38, 6.465, 5.00035E+09, -13.022, 1.993, -41.784, 13.065, -1.300, 0.04267703},
603 -5.328832253e+02, 2.398514938e+02, -4.04455 << 620 { 39, 6.492, 5.00035E+09, -13.043, 1.994, -44.609, 14.114, -1.429, 0.0479348},
604 {14, 6.000000000e+00, 7.587999300e+00, 1.499 << 621 { 40, 6.499, 5.00035E+09, -13.064, 1.994, -47.142, 15.019, -1.536, 0.0521347},
605 3.978691889e+02, -2.370975001e+02, 5.158692 << 622 { 41, 6.384, 5.00035E+09, -13.156, 1.996, -53.114, 17.052, -1.766, 0.06079426},
606 -2.340256277e+02, 9.813362251e+01, -1.52789 << 623 { 42, 6.500, 5.00035E+09, -13.176, 1.996, -54.590, 17.550, -1.822, 0.06290335},
607 {15, 6.000000000e+00, 7.587999300e+00, 1.499 << 624 { 43, 6.500, 5.00035E+09, -13.133, 1.997, -51.272, 16.423, -1.694, 0.05806108},
608 2.569833671e+02, -1.513623448e+02, 3.210087 << 625 { 44, 6.500, 5.00035E+09, -13.220, 1.996, -58.314, 18.839, -1.969, 0.0684608},
609 -1.345727293e+01, -6.291081167e+00, 3.23596 << 626 { 45, 6.500, 5.00035E+09, -13.246, 1.998, -59.674, 19.295, -2.020, 0.07037294},
610 {16, 6.000000000e+00, 7.587999300e+00, 1.499 << 627 { 46, 6.500, 5.00035E+09, -13.407, 1.999, -72.228, 23.693, -2.532, 0.09017969},
611 1.015293074e+02, -5.721639224e+01, 1.078607 << 628 { 47, 6.500, 5.00035E+09, -13.277, 1.998, -60.890, 19.647, -2.053, 0.07138694},
612 1.854818165e+02, -1.000803879e+02, 1.979815 << 629 { 48, 6.500, 5.00035E+09, -13.222, 1.998, -56.152, 18.002, -1.863, 0.06410123},
613 {17, 6.000000000e+00, 7.587999300e+00, 1.499 << 630 { 49, 6.500, 5.00035E+09, -13.199, 1.997, -56.208, 18.052, -1.872, 0.06456884},
614 -4.294163461e+01, 2.862162412e+01, -8.28597 << 631 { 50, 6.500, 5.00035E+09, -13.215, 1.998, -58.478, 18.887, -1.973, 0.06860356},
615 1.676674074e+02, -8.976414784e+01, 1.763329 << 632 { 51, 6.500, 5.00035E+09, -13.230, 1.998, -60.708, 19.676, -2.066, 0.07225841},
616 {18, 6.000000000e+00, 7.587999300e+00, 1.499 << 633 { 52, 6.500, 7.99834E+09, -13.246, 1.998, -63.341, 20.632, -2.180, 0.0767412},
617 -3.573422746e+01, 2.403066369e+01, -7.17361 << 634 { 53, 6.500, 5.00035E+09, -13.262, 1.998, -66.339, 21.716, -2.310, 0.08191981},
618 1.811925229e+02, -9.574636323e+01, 1.861940 << 635 { 54, 6.500, 7.99834E+09, -13.279, 1.998, -67.649, 22.151, -2.357, 0.08357825},
619 {19, 6.000000000e+00, 7.650499797e+00, 1.999 << 636 { 55, 6.500, 5.00035E+09, -12.951, 1.990, -45.302, 14.219, -1.423, 0.04712317},
620 1.263152069e+02, -8.738932892e+01, 2.109042 << 637 { 56, 6.425, 5.00035E+09, -12.882, 1.992, -39.825, 12.363, -1.214, 0.03931009},
621 9.183312428e+01, -5.232836676e+01, 1.072450 << 638 { 57, 6.466, 2.82488E+09, -12.903, 1.992, -38.952, 11.982, -1.160, 0.03681554},
622 {20, 6.000000000e+00, 7.650499797e+00, 1.999 << 639 { 58, 6.451, 5.00035E+09, -12.915, 1.993, -41.959, 13.118, -1.302, 0.04271291},
623 6.620218058e+02, -4.057504297e+02, 9.180787 << 640 { 59, 6.434, 5.00035E+09, -12.914, 1.993, -40.528, 12.555, -1.230, 0.03971407},
624 7.034138711e+01, -4.198325416e+01, 8.861351 << 641 { 60, 6.444, 5.00035E+09, -12.922, 1.992, -39.986, 12.329, -1.200, 0.03843737},
625 {21, 6.000000000e+00, 7.650499797e+00, 1.999 << 642 { 61, 6.414, 7.99834E+09, -12.930, 1.993, -42.756, 13.362, -1.327, 0.0436124},
626 6.766093786e+02, -4.129087029e+02, 9.305090 << 643 { 62, 6.420, 7.99834E+09, -12.938, 1.992, -42.682, 13.314, -1.319, 0.04322509},
627 1.916559096e+01, -1.807294109e+01, 4.677205 << 644 { 63, 6.416, 7.99834E+09, -12.946, 1.993, -42.399, 13.185, -1.301, 0.04243861},
628 {22, 6.000000000e+00, 7.650499797e+00, 1.999 << 645 { 64, 6.443, 7.99834E+09, -12.963, 1.993, -43.226, 13.475, -1.335, 0.04377341},
629 6.969823082e+02, -4.236620289e+02, 9.513714 << 646 { 65, 6.449, 7.99834E+09, -12.973, 1.993, -43.232, 13.456, -1.330, 0.04347536},
630 -6.501317146e+01, 2.138553133e+01, -2.25099 << 647 { 66, 6.419, 7.99834E+09, -12.966, 1.993, -42.047, 12.990, -1.270, 0.04095499},
631 {23, 6.000000000e+00, 7.650499797e+00, 1.999 << 648 { 67, 6.406, 7.99834E+09, -12.976, 1.993, -42.405, 13.106, -1.283, 0.04146024},
632 6.889749928e+02, -4.181421624e+02, 9.373529 << 649 { 68, 6.424, 7.99834E+09, -12.986, 1.993, -41.974, 12.926, -1.259, 0.040435},
633 -1.382770534e+02, 5.540647456e+01, -8.17001 << 650 { 69, 6.417, 7.99834E+09, -12.989, 1.993, -42.132, 12.967, -1.262, 0.04048908},
634 {24, 6.000000000e+00, 7.650499797e+00, 1.999 << 651 { 70, 6.405, 7.99834E+09, -13.000, 1.994, -42.582, 13.122, -1.280, 0.04119599},
635 4.365566411e+02, -2.672774427e+02, 6.001631 << 652 { 71, 6.449, 7.99834E+09, -13.015, 1.994, -42.586, 13.115, -1.278, 0.04107587},
636 -2.393534124e+02, 1.020845165e+02, -1.62474 << 653 { 72, 6.465, 7.99834E+09, -13.030, 1.994, -43.708, 13.509, -1.324, 0.04286491},
637 {25, 6.000000000e+00, 7.650499797e+00, 1.999 << 654 { 73, 6.447, 7.99834E+09, -13.048, 1.996, -44.838, 13.902, -1.369, 0.04457132},
638 6.461381990e+02, -3.918546518e+02, 8.769548 << 655 { 74, 6.452, 7.99834E+09, -13.073, 1.997, -45.545, 14.137, -1.395, 0.04553459},
639 -2.597409979e+02, 1.113332866e+02, -1.78212 << 656 { 75, 6.432, 7.99834E+09, -13.082, 1.997, -46.426, 14.431, -1.428, 0.04678218},
640 {26, 6.000000000e+00, 7.849500202e+00, 5.000 << 657 { 76, 6.439, 7.99834E+09, -13.100, 1.997, -47.513, 14.806, -1.471, 0.04842566},
641 4.261007401e+02, -2.588846763e+02, 5.764613 << 658 { 77, 6.432, 7.99834E+09, -13.110, 1.997, -48.225, 15.042, -1.497, 0.04938364},
642 -1.982896712e+02, 8.274273985e+01, -1.28407 << 659 { 78, 6.500, 7.99834E+09, -13.185, 1.997, -53.256, 16.739, -1.687, 0.05645173},
643 {27, 6.000000000e+00, 7.849500202e+00, 5.000 << 660 { 79, 6.500, 7.99834E+09, -13.200, 1.997, -53.900, 16.946, -1.709, 0.05723134},
644 4.006816638e+02, -2.439311564e+02, 5.435031 << 661 { 80, 6.500, 7.99834E+09, -13.156, 1.998, -49.801, 15.536, -1.547, 0.05103522},
645 -2.205075564e+02, 9.262919772e+01, -1.44890 << 662 { 81, 6.500, 7.99834E+09, -13.128, 1.997, -49.651, 15.512, -1.548, 0.05123203},
646 {28, 6.000000000e+00, 7.849500202e+00, 5.000 << 663 { 82, 6.500, 7.99834E+09, -13.134, 1.997, -51.021, 16.018, -1.609, 0.05364831},
647 3.967750019e+02, -2.411866801e+02, 5.364872 << 664 { 83, 6.500, 7.99834E+09, -13.148, 1.998, -52.693, 16.612, -1.679, 0.05638698},
648 -2.516823030e+02, 1.065117131e+02, -1.68053 << 665 { 84, 6.500, 7.99834E+09, -13.161, 1.998, -54.415, 17.238, -1.754, 0.05935566},
649 {29, 6.000000000e+00, 7.849500202e+00, 5.000 << 666 { 85, 6.500, 7.99834E+09, -13.175, 1.998, -56.083, 17.834, -1.824, 0.06206068},
650 2.437671888e+02, -1.499592208e+02, 3.332221 << 667 { 86, 6.500, 7.99834E+09, -13.189, 1.998, -57.860, 18.463, -1.898, 0.0649633},
651 -2.874130637e+02, 1.223381969e+02, -1.94317 << 668 { 87, 6.500, 7.99834E+09, -12.885, 1.990, -39.973, 12.164, -1.162, 0.0364598},
652 {30, 6.000000000e+00, 7.849500202e+00, 5.000 << 669 { 88, 6.417, 7.99834E+09, -12.816, 1.991, -34.591, 10.338, -0.956, 0.0287409},
653 3.914867984e+02, -2.378147085e+02, 5.284517 << 670 { 89, 6.442, 7.99834E+09, -12.831, 1.992, -36.002, 10.867, -1.021, 0.03136835},
654 -3.235063319e+02, 1.384252948e+02, -2.21184 << 671 { 90, 6.463, 7.99834E+09, -12.850, 1.993, -37.660, 11.475, -1.095, 0.03435334},
655 {31, 6.000000000e+00, 7.849500202e+00, 5.000 << 672 { 91, 6.447, 7.99834E+09, -12.852, 1.993, -37.268, 11.301, -1.071, 0.0330539},
656 4.325820127e+02, -2.614587597e+02, 5.793273 << 673 { 92, 6.439, 7.99834E+09, -12.858, 1.993, -37.695, 11.438, -1.085, 0.03376669},
657 -3.359152840e+02, 1.437507638e+02, -2.29745 << 674 { 93, 6.437, 1.00000E+10, -12.866, 1.993, -39.010, 11.927, -1.146, 0.03630848},
658 {32, 6.000000000e+00, 7.849500202e+00, 5.000 << 675 { 94, 6.432, 7.99834E+09, -12.862, 1.993, -37.192, 11.229, -1.057, 0.0325621},
659 4.388195965e+02, -2.642662297e+02, 5.834159 << 676 { 95, 6.435, 7.99834E+09, -12.869, 1.993, -37.589, 11.363, -1.072, 0.03312393},
660 -3.430730654e+02, 1.467102631e+02, -2.34316 << 677 { 96, 6.449, 1.00000E+10, -12.886, 1.993, -39.573, 12.095, -1.162, 0.03680527},
661 {33, 6.000000000e+00, 7.849500202e+00, 5.000 << 678 { 97, 6.446, 1.00000E+10, -12.892, 1.993, -40.007, 12.242, -1.178, 0.03737377},
662 3.931399547e+02, -2.363700718e+02, 5.197696 << 679 { 98, 6.421, 1.00000E+10, -12.887, 1.993, -39.509, 12.041, -1.152, 0.03629023},
663 -3.501570134e+02, 1.497141578e+02, -2.39088 << 680 { 99, 6.414, 1.00000E+10, -12.894, 1.993, -39.939, 12.183, -1.168, 0.03690464},
664 {34, 6.000000000e+00, 7.849500202e+00, 5.000 << 681 {100, 6.412, 1.00000E+10, -12.900, 1.993, -39.973, 12.180, -1.166, 0.036773}
665 3.772588127e+02, -2.256347960e+02, 4.929790 << 682 };
666 -3.481053019e+02, 1.486490112e+02, -2.37074 <<
667 {35, 6.000000000e+00, 7.849500202e+00, 5.000 <<
668 3.344685842e+02, -1.994816236e+02, 4.332267 <<
669 -3.227660675e+02, 1.370301996e+02, -2.17154 <<
670 {36, 6.000000000e+00, 7.849500202e+00, 5.000 <<
671 3.004054446e+02, -1.781334135e+02, 3.834850 <<
672 -2.980827664e+02, 1.257508661e+02, -1.97879 <<
673 {37, 6.000000000e+00, 7.849500202e+00, 5.000 <<
674 -3.687188343e+01, 1.054409719e+01, -8.51658 <<
675 -2.699384784e+02, 1.129635316e+02, -1.76144 <<
676 {38, 6.000000000e+00, 7.849500202e+00, 5.000 <<
677 1.969969064e+02, -1.286503864e+02, 3.008431 <<
678 -2.331258613e+02, 9.627987243e+01, -1.47851 <<
679 {39, 6.000000000e+00, 7.849500202e+00, 5.000 <<
680 2.891710763e+02, -1.819536752e+02, 4.158265 <<
681 -1.997404800e+02, 8.119476676e+01, -1.22342 <<
682 {40, 6.000000000e+00, 7.849500202e+00, 5.000 <<
683 3.393782172e+02, -2.103908454e+02, 4.758278 <<
684 -1.549247582e+02, 6.091403935e+01, -8.79930 <<
685 {41, 6.000000000e+00, 7.849500202e+00, 5.000 <<
686 2.748604341e+02, -1.706429616e+02, 3.843757 <<
687 -1.163607425e+02, 4.350905533e+01, -5.85930 <<
688 {42, 6.000000000e+00, 7.849500202e+00, 5.000 <<
689 3.203285955e+02, -1.966282865e+02, 4.398204 <<
690 -9.364181222e+01, 3.329814493e+01, -4.14168 <<
691 {43, 6.000000000e+00, 7.849500202e+00, 5.000 <<
692 4.184977165e+02, -2.552902161e+02, 5.707764 <<
693 -8.395646154e+01, 2.898228589e+01, -3.42235 <<
694 {44, 6.000000000e+00, 7.849500202e+00, 5.000 <<
695 3.243555305e+02, -1.978255470e+02, 4.397580 <<
696 -5.506292375e+01, 1.599310639e+01, -1.23715 <<
697 {45, 6.000000000e+00, 7.849500202e+00, 5.000 <<
698 3.037823599e+02, -1.856628295e+02, 4.128167 <<
699 -5.014186072e+01, 1.386962969e+01, -8.95080 <<
700 {46, 6.000000000e+00, 7.849500202e+00, 5.000 <<
701 3.529797051e+02, -2.101512262e+02, 4.563946 <<
702 -4.815922691e+01, 1.301508788e+01, -7.58085 <<
703 {47, 6.000000000e+00, 7.849500202e+00, 5.000 <<
704 3.074953924e+02, -1.872462583e+02, 4.149827 <<
705 -4.897188379e+01, 1.335300002e+01, -8.11005 <<
706 {48, 6.000000000e+00, 7.849500202e+00, 5.000 <<
707 4.059717166e+02, -2.462737702e+02, 5.472040 <<
708 -5.901534554e+01, 1.791385249e+01, -1.58706 <<
709 {49, 6.000000000e+00, 7.849500202e+00, 5.000 <<
710 4.369774251e+02, -2.639721687e+02, 5.849617 <<
711 -7.399698219e+01, 2.469785523e+01, -2.73788 <<
712 {50, 6.000000000e+00, 7.849500202e+00, 5.000 <<
713 4.289361021e+02, -2.585593024e+02, 5.714058 <<
714 -9.269047286e+01, 3.314422349e+01, -4.16734 <<
715 {51, 6.000000000e+00, 7.849500202e+00, 5.000 <<
716 3.866985836e+02, -2.328379698e+02, 5.128884 <<
717 -1.067869310e+02, 3.950715983e+01, -5.24332 <<
718 {52, 6.000000000e+00, 7.951499931e+00, 7.998 <<
719 3.947511198e+02, -2.363799049e+02, 5.179393 <<
720 -1.069681982e+02, 3.995521754e+01, -5.38207 <<
721 {53, 6.000000000e+00, 7.849500202e+00, 5.000 <<
722 3.694394448e+02, -2.204699428e+02, 4.806381 <<
723 -1.180749905e+02, 4.460080701e+01, -6.10521 <<
724 {54, 6.000000000e+00, 7.951499931e+00, 7.998 <<
725 3.423943987e+02, -2.041330669e+02, 4.437639 <<
726 -1.288973984e+02, 4.985324046e+01, -7.05604 <<
727 {55, 6.000000000e+00, 7.849500202e+00, 5.000 <<
728 -7.663422017e+01, 3.462700567e+01, -6.27355 <<
729 -1.318428276e+02, 5.081036112e+01, -7.15490 <<
730 {56, 6.000000000e+00, 7.849500202e+00, 5.000 <<
731 1.084179205e+02, -7.602229206e+01, 1.843754 <<
732 -1.346311376e+02, 5.207427468e+01, -7.36983 <<
733 {57, 6.000000000e+00, 7.725500002e+00, 2.824 <<
734 2.995898890e+02, -1.889477671e+02, 4.336642 <<
735 5.503972208e+00, -1.227641064e+01, 3.699182 <<
736 {58, 6.000000000e+00, 7.849500202e+00, 5.000 <<
737 1.709135500e+02, -1.120124681e+02, 2.615893 <<
738 -1.375860132e+02, 5.337811974e+01, -7.58678 <<
739 {59, 6.000000000e+00, 7.849500202e+00, 5.000 <<
740 1.214691988e+02, -8.336119630e+01, 1.996468 <<
741 -1.631005912e+02, 6.472051894e+01, -9.47609 <<
742 {60, 6.000000000e+00, 7.849500202e+00, 5.000 <<
743 1.302719596e+02, -8.835087414e+01, 2.101971 <<
744 -1.692901279e+02, 6.742727614e+01, -9.92066 <<
745 {61, 6.000000000e+00, 7.951499931e+00, 7.998 <<
746 1.127680235e+02, -7.782238836e+01, 1.865126 <<
747 -2.059821608e+02, 8.384774285e+01, -1.26734 <<
748 {62, 6.000000000e+00, 7.951499931e+00, 7.998 <<
749 1.203145109e+02, -8.212556537e+01, 1.956606 <<
750 -2.158058793e+02, 8.810144391e+01, -1.33638 <<
751 {63, 6.000000000e+00, 7.951499931e+00, 7.998 <<
752 1.212159597e+02, -8.256559477e+01, 1.964122 <<
753 -2.278531434e+02, 9.336519465e+01, -1.42258 <<
754 {64, 6.000000000e+00, 7.951499931e+00, 7.998 <<
755 1.689382403e+02, -1.099987696e+02, 2.551961 <<
756 -2.282716670e+02, 9.348611199e+01, -1.42358 <<
757 {65, 6.000000000e+00, 7.951499931e+00, 7.998 <<
758 1.724155378e+02, -1.120798437e+02, 2.598264 <<
759 -2.322687147e+02, 9.517466656e+01, -1.45033 <<
760 {66, 6.000000000e+00, 7.951499931e+00, 7.998 <<
761 1.286079419e+02, -8.646296410e+01, 2.039801 <<
762 -2.420048480e+02, 9.935663043e+01, -1.51765 <<
763 {67, 6.000000000e+00, 7.951499931e+00, 7.998 <<
764 1.182799697e+02, -8.043389241e+01, 1.908027 <<
765 -2.464462609e+02, 1.012059056e+02, -1.54646 <<
766 {68, 6.000000000e+00, 7.951499931e+00, 7.998 <<
767 1.150510247e+02, -7.859576077e+01, 1.868688 <<
768 -2.457555063e+02, 1.007538481e+02, -1.53669 <<
769 {69, 6.000000000e+00, 7.951499931e+00, 7.998 <<
770 1.266280406e+02, -8.514491730e+01, 2.007089 <<
771 -2.492442707e+02, 1.021615320e+02, -1.55787 <<
772 {70, 6.000000000e+00, 7.951499931e+00, 7.998 <<
773 1.224253568e+02, -8.281395858e+01, 1.958609 <<
774 -2.488808342e+02, 1.018569466e+02, -1.55060 <<
775 {71, 6.000000000e+00, 7.951499931e+00, 7.998 <<
776 1.862181262e+02, -1.199038630e+02, 2.763107 <<
777 -2.403102476e+02, 9.796272016e+01, -1.48452 <<
778 {72, 6.000000000e+00, 7.951499931e+00, 7.998 <<
779 2.297759959e+02, -1.448485621e+02, 3.295877 <<
780 -2.282155654e+02, 9.249921555e+01, -1.39226 <<
781 {73, 6.000000000e+00, 7.951499931e+00, 7.998 <<
782 2.646909006e+02, -1.647716545e+02, 3.719903 <<
783 -2.165150972e+02, 8.722660467e+01, -1.30341 <<
784 {74, 6.000000000e+00, 7.951499931e+00, 7.998 <<
785 2.251239174e+02, -1.414731209e+02, 3.206048 <<
786 -2.070173544e+02, 8.296725365e+01, -1.23198 <<
787 {75, 6.000000000e+00, 7.951499931e+00, 7.998 <<
788 2.627532736e+02, -1.629008146e+02, 3.661592 <<
789 -1.945762063e+02, 7.740995255e+01, -1.13912 <<
790 {76, 6.000000000e+00, 7.951499931e+00, 7.998 <<
791 2.644549626e+02, -1.637369900e+02, 3.675734 <<
792 -1.725967865e+02, 6.755389456e+01, -9.73763 <<
793 {77, 6.000000000e+00, 7.951499931e+00, 7.998 <<
794 2.677629012e+02, -1.650589135e+02, 3.690999 <<
795 -1.584140848e+02, 6.122430396e+01, -8.68087 <<
796 {78, 6.000000000e+00, 7.951499931e+00, 7.998 <<
797 2.420702029e+02, -1.484461630e+02, 3.292288 <<
798 -1.319886050e+02, 4.940494114e+01, -6.70274 <<
799 {79, 6.000000000e+00, 7.951499931e+00, 7.998 <<
800 2.346714957e+02, -1.439356552e+02, 3.189416 <<
801 -1.130109430e+02, 4.093029258e+01, -5.28674 <<
802 {80, 6.000000000e+00, 7.951499931e+00, 7.998 <<
803 2.747370538e+02, -1.689673404e+02, 3.771696 <<
804 -9.001823908e+01, 3.066094857e+01, -3.57045 <<
805 {81, 6.000000000e+00, 7.951499931e+00, 7.998 <<
806 3.142563781e+02, -1.916613838e+02, 4.259167 <<
807 -7.642731867e+01, 2.462410146e+01, -2.56697 <<
808 {82, 6.000000000e+00, 7.951499931e+00, 7.998 <<
809 3.509258060e+02, -2.125470710e+02, 4.702461 <<
810 -5.173355302e+01, 1.362015056e+01, -7.32128 <<
811 {83, 6.000000000e+00, 7.951499931e+00, 7.998 <<
812 3.399729483e+02, -2.056319770e+02, 4.539614 <<
813 -4.131443229e+01, 8.986236911e+00, 3.924628 <<
814 {84, 6.000000000e+00, 7.951499931e+00, 7.998 <<
815 3.640602841e+02, -2.190164327e+02, 4.815603 <<
816 -3.256862965e+01, 5.115606198e+00, 6.800853 <<
817 {85, 6.000000000e+00, 7.951499931e+00, 7.998 <<
818 3.766488275e+02, -2.257321142e+02, 4.947300 <<
819 -2.300947210e+01, 8.615223509e-01, 1.388425 <<
820 {86, 6.000000000e+00, 7.951499931e+00, 7.998 <<
821 3.443622947e+02, -2.064342780e+02, 4.516044 <<
822 -5.399039282e+00, -7.002814559e+00, 2.70251 <<
823 {87, 6.000000000e+00, 7.951499931e+00, 7.998 <<
824 -3.706791591e+01, 1.118013187e+01, -1.05772 <<
825 -3.451314336e+00, -7.779254134e+00, 2.81626 <<
826 {88, 6.000000000e+00, 7.951499931e+00, 7.998 <<
827 6.125934670e+01, -4.855548659e+01, 1.248551 <<
828 -6.021643455e+00, -6.580234329e+00, 2.60744 <<
829 {89, 6.000000000e+00, 7.951499931e+00, 7.998 <<
830 1.350863292e+02, -9.126618691e+01, 2.169932 <<
831 1.937135880e+01, -1.787129621e+01, 4.485878 <<
832 {90, 6.000000000e+00, 7.951499931e+00, 7.998 <<
833 1.784388998e+02, -1.161623817e+02, 2.702376 <<
834 2.216057166e+01, -1.904990091e+01, 4.671627 <<
835 {91, 6.000000000e+00, 7.951499931e+00, 7.998 <<
836 1.368355213e+02, -9.179790820e+01, 2.169910 <<
837 4.516580666e+00, -1.118102949e+01, 3.357662 <<
838 {92, 6.000000000e+00, 7.951499931e+00, 7.998 <<
839 1.427130850e+02, -9.499714618e+01, 2.234475 <<
840 1.341991149e+01, -1.518503354e+01, 4.030838 <<
841 {93, 6.000000000e+00, 8.000000000e+00, 1.000 <<
842 2.341801100e+01, -2.506119713e+01, 7.023029 <<
843 -3.575331738e+01, 7.276302226e+00, 1.906771 <<
844 {94, 6.000000000e+00, 7.951499931e+00, 7.998 <<
845 1.287618322e+02, -8.721780968e+01, 2.073255 <<
846 -2.307262580e+01, 1.113132278e+00, 1.305250 <<
847 {95, 6.000000000e+00, 7.951499931e+00, 7.998 <<
848 1.334821220e+02, -8.985337775e+01, 2.127928 <<
849 -3.518662723e+01, 6.514543434e+00, 4.030862 <<
850 {96, 6.000000000e+00, 8.000000000e+00, 1.000 <<
851 4.545581472e+01, -3.771304300e+01, 9.729129 <<
852 -4.973805034e+01, 1.342335334e+01, -8.22113 <<
853 {97, 6.000000000e+00, 8.000000000e+00, 1.000 <<
854 4.689042092e+01, -3.843347264e+01, 9.859294 <<
855 -4.657434145e+01, 1.204637835e+01, -5.98244 <<
856 {98, 6.000000000e+00, 8.000000000e+00, 1.000 <<
857 1.337584189e+01, -1.907284620e+01, 5.691614 <<
858 -5.573362773e+01, 1.615667599e+01, -1.28896 <<
859 {99, 6.000000000e+00, 8.000000000e+00, 1.000 <<
860 1.376201293e+01, -1.919251815e+01, 5.693799 <<
861 -4.914211254e+01, 1.314247998e+01, -7.73933 <<
862 {100, 6.000000000e+00, 8.000000000e+00, 1.00 <<
863 1.277081775e+01, -1.854047224e+01, 5.534680 <<
864 -5.074293980e+01, 1.383260974e+01, -8.85890 <<
865 //....oooOO0OOooo........oooOO0OOooo........oo <<
866 683