ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4AdjointCSMatrix.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4AdjointCSMatrix.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 
28 #include <iomanip>
29 #include <fstream>
30 
31 #include "G4AdjointCSMatrix.hh"
32 #include "G4SystemOfUnits.hh"
33 #include "G4AdjointInterpolator.hh"
34 
36 //
37 G4AdjointCSMatrix::G4AdjointCSMatrix(G4bool aBool){
38  theLogPrimEnergyVector.clear();
39  theLogCrossSectionVector.clear();
40  theLogSecondEnergyMatrix.clear();
41  theLogProbMatrix.clear();
42  theLogProbMatrixIndex.clear();
43  log0Vector.clear();
44  nb_of_PrimEnergy=0;
45  is_scat_proj_to_proj_case =aBool;
46  dlog =0;
47 }
49 //
50 G4AdjointCSMatrix::~G4AdjointCSMatrix(){
51  theLogPrimEnergyVector.clear();
52  theLogCrossSectionVector.clear();
53  theLogSecondEnergyMatrix.clear();
54  theLogProbMatrix.clear();
55 }
57 //
58 void G4AdjointCSMatrix::Clear()
59 {
60  theLogPrimEnergyVector.clear();
61  theLogCrossSectionVector.clear();
62  theLogSecondEnergyMatrix.clear();
63  theLogProbMatrix.clear();
64  theLogProbMatrixIndex.clear();
65  log0Vector.clear();
66  nb_of_PrimEnergy=0;
67 }
69 //
70  void G4AdjointCSMatrix::AddData(G4double aLogPrimEnergy,G4double aLogCS, std::vector< double>* aLogSecondEnergyVector,
71  std::vector< double>* aLogProbVector,size_t n_pro_decade){
72 
73  G4AdjointInterpolator* theInterpolator=G4AdjointInterpolator::GetInstance();
74 
75  //At this time we consider that the energy is increasing monotically
76  theLogPrimEnergyVector.push_back(aLogPrimEnergy);
77  theLogCrossSectionVector.push_back(aLogCS);
78  theLogSecondEnergyMatrix.push_back(aLogSecondEnergyVector);
79  theLogProbMatrix.push_back(aLogProbVector);
80 
81  std::vector< size_t>* aLogProbVectorIndex = 0;
82  dlog =0;
83 
84  if (n_pro_decade > 0 && aLogProbVector->size()>0) {
85  aLogProbVectorIndex = new std::vector< size_t>();
86  dlog=std::log(10.)/n_pro_decade;
87  G4double log_val = int(std::min((*aLogProbVector)[0],aLogProbVector->back())/dlog)*dlog;
88  log0Vector.push_back(log_val);
89 
90  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
91  while(log_val<0.) {
92  aLogProbVectorIndex->push_back(theInterpolator->FindPosition(log_val,(*aLogProbVector)));
93  log_val+=dlog;
94  }
95  }
96  else {
97  log0Vector.push_back(0.);
98  }
99  theLogProbMatrixIndex.push_back(aLogProbVectorIndex);
100 
101 
102  nb_of_PrimEnergy++;
103 
104 
105 }
107 //
108 G4bool G4AdjointCSMatrix::GetData(unsigned int i, G4double& aLogPrimEnergy,G4double& aLogCS,G4double& log0, std::vector< double>*& aLogSecondEnergyVector,
109  std::vector< double>*& aLogProbVector, std::vector< size_t>*& aLogProbVectorIndex)
110 { if (i>= nb_of_PrimEnergy) return false;
111  //G4cout<<"Test Get Data "<<G4endl;
112  aLogPrimEnergy = theLogPrimEnergyVector[i];
113  aLogCS = theLogCrossSectionVector[i];
114  aLogSecondEnergyVector = theLogSecondEnergyMatrix[i];
115  aLogProbVector = theLogProbMatrix[i];
116  aLogProbVectorIndex = theLogProbMatrixIndex[i];
117  log0=log0Vector[i];
118  return true;
119 
120 }
122 //
123 void G4AdjointCSMatrix::Write(G4String file_name)
124 { std::fstream FileOutput(file_name, std::ios::out);
125  FileOutput<<std::setiosflags(std::ios::scientific);
126  FileOutput<<std::setprecision(6);
127  FileOutput<<theLogPrimEnergyVector.size()<<G4endl;
128  for (size_t i=0;i<theLogPrimEnergyVector.size();i++){
129  FileOutput<<std::exp(theLogPrimEnergyVector[i])/MeV<<'\t'<<std::exp(theLogCrossSectionVector[i])<<G4endl;
130  size_t j1=0;
131  FileOutput<<theLogSecondEnergyMatrix[i]->size()<<G4endl;
132  for (size_t j=0;j<theLogSecondEnergyMatrix[i]->size();j++){
133  FileOutput<<std::exp((*theLogSecondEnergyMatrix[i])[j]);
134  j1++;
135  if (j1<10) FileOutput<<'\t';
136  else {
137  FileOutput<<G4endl;
138  j1=0;
139  }
140  }
141  if (j1>0) FileOutput<<G4endl;
142  j1=0;
143  FileOutput<<theLogProbMatrix[i]->size()<<G4endl;
144  for (size_t j=0;j<theLogProbMatrix[i]->size();j++){
145  FileOutput<<std::exp((*theLogProbMatrix[i])[j]);
146  j1++;
147  if (j1<10) FileOutput<<'\t';
148  else {
149  FileOutput<<G4endl;
150  j1=0;
151  }
152  }
153  if (j1>0) FileOutput<<G4endl;
154 
155 
156  }
157 
158 }
160 //
161 void G4AdjointCSMatrix::Read(G4String file_name)
162 { std::fstream FileOutput(file_name, std::ios::in);
163  size_t n1,n2;
164 
165 
166  theLogPrimEnergyVector.clear();
167  theLogCrossSectionVector.clear();
168  theLogSecondEnergyMatrix.clear();
169  theLogProbMatrix.clear();
170  FileOutput>>n1;
171  for (size_t i=0; i<n1;i++){
172  G4double E,CS;
173  FileOutput>>E>>CS;
174  theLogPrimEnergyVector.push_back(E);
175  theLogCrossSectionVector.push_back(CS);
176  FileOutput>>n2;
177  theLogSecondEnergyMatrix.push_back(new std::vector<G4double>());
178  theLogProbMatrix.push_back(new std::vector<G4double>());
179 
180  for (size_t j=0; j<n2;j++){
181  G4double E1;
182  FileOutput>>E1;
183  theLogSecondEnergyMatrix[i]->push_back(E1);
184  }
185  FileOutput>>n2;
186  for (size_t j=0; j<n2;j++){
187  G4double prob;
188  FileOutput>>prob;
189  theLogProbMatrix[i]->push_back(prob);
190  }
191 
192 
193 
194  }
195 
196 
197 
198 
199 }