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