GENIEGenerator
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
PathLengthList.cxx
Go to the documentation of this file.
1 //____________________________________________________________________________
2 /*
3  Copyright (c) 2003-2024, The GENIE Collaboration
4  For the full text of the license visit http://copyright.genie-mc.org
5 
6  Costas Andreopoulos <c.andreopoulos \at cern.ch>
7  University of Liverpool
8 */
9 //____________________________________________________________________________
10 
11 #include <fstream>
12 #include <iomanip>
13 
14 #include "libxml/parser.h"
15 #include "libxml/xmlmemory.h"
16 
17 #include <TLorentzVector.h>
18 
27 
28 
29 using std::ofstream;
30 using std::setw;
31 using std::setfill;
32 using std::endl;
33 
34 using namespace genie;
35 
36 //____________________________________________________________________________
37 namespace genie {
38  ostream & operator << (ostream & stream, const PathLengthList & list)
39  {
40  list.Print(stream);
41  return stream;
42  }
43 }
44 //___________________________________________________________________________
46 map<int, double>()
47 {
48 
49 }
50 //___________________________________________________________________________
52 map<int, double>()
53 {
54  PDGCodeList::const_iterator pdg_iter;
55 
56  for(pdg_iter = pdg_list.begin(); pdg_iter != pdg_list.end(); ++pdg_iter) {
57  int pdgc = *pdg_iter;
58  this->insert( map<int, double>::value_type(pdgc, 0.) );
59  }
60 }
61 //___________________________________________________________________________
63 map<int, double>()
64 {
65  this->Copy(plist);
66 }
67 //__________________________________________________________________________
68 PathLengthList::PathLengthList(const map<int,double> & plist) :
69 map<int, double>()
70 {
71  map<int,double>::const_iterator iter;
72 
73  for(iter = plist.begin(); iter != plist.end(); ++iter) {
74  int pdgc = iter->first;
75  double pl = iter->second;
76  this->insert( map<int, double>::value_type(pdgc, pl) );
77  }
78 }
79 //__________________________________________________________________________
81 {
82 
83 }
84 //___________________________________________________________________________
85 void PathLengthList::AddPathLength(int pdgc, double pl)
86 {
87 // Adds pl to the total path length for material with code = pdgc
88 
89  if (this->count(pdgc) == 1) { (*this)[pdgc] += pl; }
90  else {
91  LOG("PathL", pWARN)
92  << "No material with PDG code = " << pdgc << " in path length list";
93  }
94 }
95 //___________________________________________________________________________
96 void PathLengthList::SetPathLength(int pdgc, double pl)
97 {
98 // Sets the total path length for material with code = pdgc to be pl
99 
100  if (this->count(pdgc) == 1) { (*this)[pdgc] = pl; }
101  else {
102  LOG("PathL", pWARN)
103  << "No material with PDG code = " << pdgc << " in path length list";
104  }
105 }
106 //___________________________________________________________________________
107 void PathLengthList::ScalePathLength(int pdgc, double scale)
108 {
109 // Scales pl for material with code = pdgc with the input scale factor
110 
111  if (this->count(pdgc) == 1) {
112  double pl = (*this)[pdgc];
113  pl *= scale;
114  (*this)[pdgc] = pl;
115  } else {
116  LOG("PathL", pWARN)
117  << "No material with PDG code = " << pdgc << " in path length list";
118  }
119 }
120 //___________________________________________________________________________
121 double PathLengthList::PathLength(int pdgc) const
122 {
123 // Gets the total path length for material with code = pdgc to be pl
124 
125  if ( this->count(pdgc) == 1 ) {
126  map<int, double>::const_iterator pl_iter = this->find(pdgc);
127  return pl_iter->second;
128  } else {
129  LOG("PathL", pWARN)
130  << "No material with PDG code = " << pdgc << " in path length list";
131  }
132  return 0;
133 }
134 //___________________________________________________________________________
136 {
137  PathLengthList::const_iterator pl_iter;
138 
139  for(pl_iter = this->begin(); pl_iter != this->end(); ++pl_iter) {
140  int pdgc = pl_iter->first;
141  (*this)[pdgc] = 0.;
142  }
143 }
144 //___________________________________________________________________________
146 {
147  bool allzero = true;
148 
149  PathLengthList::const_iterator pl_iter;
150 
151  for(pl_iter = this->begin(); pl_iter != this->end(); ++pl_iter) {
152  double pl = pl_iter->second;
153  allzero = allzero && (utils::math::AreEqual(pl,0.));
154  }
155  return allzero;
156 }
157 //___________________________________________________________________________
159 {
160  this->clear();
161  PathLengthList::const_iterator pl_iter;
162  for(pl_iter = plist.begin(); pl_iter != plist.end(); ++pl_iter) {
163  int pdgc = pl_iter->first;
164  double pl = pl_iter->second;
165  this->insert( map<int, double>::value_type(pdgc, pl) );
166  }
167 }
168 //___________________________________________________________________________
169 void PathLengthList::Print(ostream & stream) const
170 {
171  stream << "\n[-]" << endl;
172 
173  PDGLibrary * pdglib = PDGLibrary::Instance();
174 
175  PathLengthList::const_iterator pl_iter;
176  size_t nc = this->size();
177 
178  for(pl_iter = this->begin(); pl_iter != this->end(); ++pl_iter) {
179 
180  int pdgc = pl_iter->first;
181  double pl = pl_iter->second; // path length
182 
183  TParticlePDG * p = pdglib->Find(pdgc);
184 
185  if(!p) {
186  stream << " |---o ** ERR: no particle with PDG code: " << pdgc;
187  } else {
188  string name = p->GetName();
189  stream << " |---o code: " << pdgc << " [" << setfill(' ')
190  << setw(5) << name << "] " << "-----> path-length = " << pl;
191  }
192  if( (--nc) > 0) stream << endl;
193  }
194 }
195 //___________________________________________________________________________
197 {
198  this->clear();
199  PDGLibrary * pdglib = PDGLibrary::Instance();
200 
201  LOG("PathL", pINFO)
202  << "Loading PathLengthList from XML file: " << filename;
203 
204  xmlDocPtr xml_doc = xmlParseFile(filename.c_str() );
205 
206  if(xml_doc==NULL) {
207  LOG("PathL", pERROR)
208  << "XML file could not be parsed! [filename: " << filename << "]";
209  return kXmlNotParsed;
210  }
211 
212  xmlNodePtr xmlCur = xmlDocGetRootElement(xml_doc);
213 
214  if(xmlCur==NULL) {
215  LOG("PathL", pERROR)
216  << "XML doc. has null root element! [filename: " << filename << "]";
217  return kXmlEmpty;
218  }
219 
220  if( xmlStrcmp(xmlCur->name, (const xmlChar *) "path_length_list") ) {
221  LOG("PathL", pERROR)
222  << "XML doc. has invalid root element! [filename: " << filename << "]";
223  return kXmlInvalidRoot;
224  }
225 
226  LOG("PathL", pINFO) << "XML file was successfully parsed";
227 
228  xmlCur = xmlCur->xmlChildrenNode; // <path_length>'s
229 
230  // loop over all xml tree nodes that are children of the root node
231  while (xmlCur != NULL) {
232 
233  // enter everytime you find a <path_length> tag
234  if( (!xmlStrcmp(xmlCur->name, (const xmlChar *) "path_length")) ) {
235 
236  xmlNodePtr xmlPlVal = xmlCur->xmlChildrenNode;
237 
238  string spdgc = utils::str::TrimSpaces(
239  utils::xml::GetAttribute(xmlCur, "pdgc"));
240 
241  string spl = utils::xml::TrimSpaces(
242  xmlNodeListGetString(xml_doc, xmlPlVal, 1));
243 
244  LOG("PathL", pDEBUG) << "pdgc = " << spdgc << " --> pl = " << spl;
245 
246  int pdgc = atoi( spdgc.c_str() );
247  double pl = atof( spl.c_str() );
248 
249  TParticlePDG * p = pdglib->Find(pdgc);
250  if(!p) {
251  LOG("PathL", pERROR)
252  << "No particle with pdgc " << pdgc
253  << " found. Will not load its path length";
254  } else
255  this->insert( map<int, double>::value_type(pdgc, pl) );
256 
257  xmlFree(xmlPlVal);
258  }
259  xmlCur = xmlCur->next;
260  } // [end of] loop over tags within root elements
261 
262  xmlFree(xmlCur);
263  return kXmlOK;
264 }
265 //___________________________________________________________________________
266 void PathLengthList::SaveAsXml(string filename) const
267 {
268 //! Save path length list to XML file
269 
270  LOG("PathL", pINFO)
271  << "Saving PathLengthList as XML in file: " << filename;
272 
273  PDGLibrary * pdglib = PDGLibrary::Instance();
274 
275  ofstream outxml(filename.c_str());
276  if(!outxml.is_open()) {
277  LOG("PathL", pERROR) << "Couldn't create file = " << filename;
278  return;
279  }
280  outxml << "<?xml version=\"1.0\" encoding=\"ISO-8859-1\"?>";
281  outxml << endl << endl;
282  outxml << "<!-- generated by PathLengthList::SaveAsXml() -->";
283  outxml << endl << endl;
284 
285  outxml << "<path_length_list>" << endl << endl;
286 
287  PathLengthList::const_iterator pl_iter;
288 
289  for(pl_iter = this->begin(); pl_iter != this->end(); ++pl_iter) {
290 
291  int pdgc = pl_iter->first;
292  double pl = pl_iter->second; // path length
293 
294  TParticlePDG * p = pdglib->Find(pdgc);
295 
296  outxml << " <path_length pdgc=\"" << pdgc << "\"> "
297  << setfill(' ') << setw(10) << pl << " </path_length>";
298  if ( p ) outxml << " <!-- [" << setfill(' ')
299  << setw(5) << p->GetName() << "] -->";
300  outxml << endl;
301  }
302  outxml << endl << "</path_length_list>";
303  outxml << endl;
304 
305  outxml.close();
306 
307 }
308 //___________________________________________________________________________
310 {
311  this->Copy(list);
312  return (*this);
313 }
314 //___________________________________________________________________________
void ScalePathLength(int pdgc, double scale)
double PathLength(int pdgc) const
bool AreAllZero(void) const
#define pERROR
Definition: Messenger.h:59
PathLengthList & operator=(const PathLengthList &list)
bool AreEqual(double x1, double x2)
Definition: MathUtils.cxx:236
string TrimSpaces(xmlChar *xmls)
void SetPathLength(int pdgc, double pl)
A list of PDG codes.
Definition: PDGCodeList.h:32
Object to be filled with the neutrino path-length, for all detector geometry materials, when starting from a position x and travelling along the direction of the neutrino 4-momentum.
void SaveAsXml(string filename) const
#define LOG(stream, priority)
A macro that returns the requested log4cpp::Category appending a string (using the FILE...
Definition: Messenger.h:96
XmlParserStatus_t LoadFromXml(string filename)
#define pINFO
Definition: Messenger.h:62
#define pWARN
Definition: Messenger.h:60
string TrimSpaces(string input)
Definition: StringUtils.cxx:18
void Print(ostream &stream) const
static PDGLibrary * Instance(void)
Definition: PDGLibrary.cxx:68
Singleton class to load &amp; serve a TDatabasePDG.
Definition: PDGLibrary.h:35
ostream & operator<<(ostream &stream, const AlgConfigPool &config_pool)
vector< vector< double > > clear
void AddPathLength(int pdgc, double pl)
TParticlePDG * Find(int pdgc, bool must_exist=true)
Definition: PDGLibrary.cxx:86
enum genie::EXmlParseStatus XmlParserStatus_t
string GetAttribute(xmlNodePtr xml_cur, string attr_name)
void Copy(const PathLengthList &plist)
#define pDEBUG
Definition: Messenger.h:63