source: snovis/trunk/source/G4Lab/cxx/PhysicsTableAccessor.cxx@ 233

Last change on this file since 233 was 233, checked in by barrand, 19 years ago
  • Property svn:eol-style set to native
File size: 14.1 KB
Line 
1// this :
2#include <G4Lab/PhysicsTableAccessor.h>
3
4// Inventor :
5#include <Inventor/SbString.h>
6#include <Inventor/SbName.h>
7#include <Inventor/nodes/SoSeparator.h>
8
9#ifdef WIN32
10#undef pascal // Clash between windef.h and Geant4/SystemOfUnits.hh
11#endif
12
13// Geant4 :
14#include <G4ParticleTable.hh>
15#include <G4PhysicsTable.hh>
16#include <G4ParticleDefinition.hh>
17#include <G4Element.hh>
18#include <G4VProcess.hh>
19#include <G4ProcessVector.hh>
20#include <G4ProcessManager.hh>
21
22// Lib :
23#include <Slash/Core/ISession.h>
24#include <Slash/Data/IIterator.h>
25#include <Lib/Out.h>
26#include <Lib/Value.h>
27#include <Lib/smanip.h>
28#include <Lib/fmanip.h>
29#include <Lib/dirmanip.h>
30#include <Lib/Debug.h>
31
32// AIDA :
33#include <AIDA/IAnalysisFactory.h>
34#include <AIDA/ITreeFactory.h>
35#include <AIDA/ITree.h>
36#include <AIDA/IPlotterFactory.h>
37#include <AIDA/IPlotter.h>
38#include <AIDA/IPlotterRegion.h>
39#include <AIDA/IFunctionFactory.h>
40#include <AIDA/IFunction.h>
41
42namespace G4Lab {
43class PhysicsTable {
44public:
45 PhysicsTable(const std::string& aTableName,
46 const std::string& aProcessName,
47 const std::string& aParticleName,
48 G4PhysicsTable* aTable,
49 G4VProcess* aProcess,
50 G4ParticleDefinition* aParticle)
51 :fTableName(aTableName)
52 ,fProcessName(aProcessName)
53 ,fParticleName(aParticleName)
54 ,fTable(aTable)
55 ,fProcess(aProcess)
56 ,fParticle(aParticle){}
57 virtual ~PhysicsTable() {
58 delete fTable;
59 }
60public:
61 std::string fTableName;
62 std::string fProcessName;
63 std::string fParticleName;
64 G4PhysicsTable* fTable;
65 G4VProcess* fProcess;
66 G4ParticleDefinition* fParticle;
67 Lib::Debug fDebug;
68};
69}
70
71//////////////////////////////////////////////////////////////////////////////
72G4Lab::PhysicsTableAccessor::PhysicsTableAccessor(
73 Slash::Core::ISession& aSession
74,AIDA::IAnalysisFactory* aAIDA
75)
76:Lib::BaseAccessor(aSession.printer())
77,fSession(aSession)
78,fType("PhysicsTable")
79,fAIDA(aAIDA)
80//////////////////////////////////////////////////////////////////////////////
81//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
82{
83 addProperty("table",Lib::Property::STRING);
84 addProperty("process",Lib::Property::STRING);
85 addProperty("particle",Lib::Property::STRING);
86 addProperty("number",Lib::Property::INTEGER);
87 addProperty("name",Lib::Property::STRING);
88}
89//////////////////////////////////////////////////////////////////////////////
90G4Lab::PhysicsTableAccessor::~PhysicsTableAccessor(
91)
92//////////////////////////////////////////////////////////////////////////////
93//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
94{
95 for(unsigned int index=0;index<fTables.size();index++) delete fTables[index];
96}
97//////////////////////////////////////////////////////////////////////////////
98std::string G4Lab::PhysicsTableAccessor::name(
99) const
100//////////////////////////////////////////////////////////////////////////////
101//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
102{
103 return fType;
104}
105namespace G4Lab {
106 class PhysicsTableIterator : public virtual Slash::Data::IIterator {
107 public: //Slash::Data::IIterator
108 virtual Slash::Data::IAccessor::Data object() {
109 if(fIterator==fVector.end()) return 0;
110 return *fIterator;
111 }
112 virtual void next() { ++fIterator; }
113 virtual void* tag() { return 0;}
114 public:
115 PhysicsTableIterator(std::vector<G4Lab::PhysicsTable*>& aVector)
116 :fVector(aVector) {
117 fIterator = fVector.begin();
118 }
119 virtual ~PhysicsTableIterator() {}
120 private:
121 std::vector<G4Lab::PhysicsTable*>& fVector;
122 std::vector<G4Lab::PhysicsTable*>::iterator fIterator;
123 };
124}
125//////////////////////////////////////////////////////////////////////////////
126Slash::Data::IIterator* G4Lab::PhysicsTableAccessor::iterator(
127)
128//////////////////////////////////////////////////////////////////////////////
129//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
130{
131 if(!fTables.size() ) {
132 G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
133 if(particleTable) {
134 G4ParticleTable::G4PTblDicIterator* particleTableIterator =
135 particleTable->GetIterator();
136 if(particleTableIterator) {
137 particleTableIterator->reset();
138 while( (*particleTableIterator)() ){
139 G4ParticleDefinition* particle = particleTableIterator->value();
140 if(particle) {
141 G4String particleName = particle->GetParticleName();
142 //Lib::Out out(fSession.printer());
143 //out << "debug : Particle : " << particleName << Lib::endl;
144 // Processes for this particle :
145 G4ProcessManager* processManager = particle->GetProcessManager();
146 if(processManager) {
147 G4ProcessVector* processList = processManager->GetProcessList();
148 if(processList) {
149 int number = processList->entries();
150 for (int index=0;index<number;index++){
151 G4VProcess* process = (*processList)(index);
152 std::vector<std::string> names;
153 std::vector<G4PhysicsTable*> tables;
154 if(findProcessTables(particle,
155 process,
156 names,
157 tables)) {
158 for(unsigned int i=0;i<tables.size();i++) {
159 std::vector<std::string> words;
160 Lib::smanip::words(names[i],".",words);
161 if(words.size()==3) {
162 fTables.push_back(new PhysicsTable(words[0],
163 words[1],
164 words[2],
165 tables[i],
166 process,
167 particle));
168 }
169 }
170 }
171 }
172 }
173 }
174 }
175 }
176 }
177 }
178 }
179
180 return new PhysicsTableIterator(fTables);
181}
182//////////////////////////////////////////////////////////////////////////////
183Slash::Core::IValue* G4Lab::PhysicsTableAccessor::findValue(
184 Slash::Data::IAccessor::Data aData
185,const std::string& aName
186,void*
187)
188//////////////////////////////////////////////////////////////////////////////
189//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
190{
191 PhysicsTable* obj = (PhysicsTable*)aData;
192 if(aName=="table") {
193 return new Lib::Value(obj->fTableName);
194 } else if(aName=="process") {
195 return new Lib::Value(obj->fProcessName);
196 } else if(aName=="particle") {
197 return new Lib::Value(obj->fParticleName);
198 } else if(aName=="number") {
199 return new Lib::Value((int)obj->fTable->size());
200 } else if(aName=="name") {
201 std::string name =
202 obj->fTableName + "." + obj->fProcessName + "." + obj->fParticleName;
203 return new Lib::Value(name);
204 } else {
205 return new Lib::Value();
206 }
207}
208//////////////////////////////////////////////////////////////////////////////
209void G4Lab::PhysicsTableAccessor::beginVisualize(
210)
211//////////////////////////////////////////////////////////////////////////////
212//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
213{
214 if(!fAIDA) {
215 Lib::Out out(fSession.printer());
216 out << "G4Lab::PhysicsTableAccessor::beginVisualize :"
217 << " AIDA not found. Can't then find an AIDA plotter."
218 << Lib::endl;
219 }
220}
221//////////////////////////////////////////////////////////////////////////////
222void G4Lab::PhysicsTableAccessor::visualize(
223 Slash::Data::IAccessor::Data aData
224,void*
225)
226//////////////////////////////////////////////////////////////////////////////
227//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
228{
229 PhysicsTable* obj = (PhysicsTable*)aData;
230
231 std::string value;
232
233 bool superpose = true;
234 if(fSession.parameterValue("modeling.superpose",value))
235 if(!Lib::smanip::tobool(value,superpose)) superpose = true;
236
237 if(!fAIDA) return;
238
239 AIDA::ITreeFactory* tf = fAIDA->createTreeFactory();
240 if(!tf) return;
241 AIDA::ITree* memoryTree = tf->create();
242 delete tf;
243 if(!memoryTree) return;
244
245 AIDA::IFunctionFactory* functionFactory =
246 fAIDA->createFunctionFactory(*memoryTree);
247 if(!functionFactory) {
248 Lib::Out out(fSession.printer());
249 out << "Can't create an function factory." << Lib::endl;
250 return;
251 }
252
253 AIDA::IPlotterFactory* plotterFactory = fAIDA->createPlotterFactory();
254 if(!plotterFactory) {
255 Lib::Out out(fSession.printer());
256 out << "Can't create a plotter factory." << Lib::endl;
257 delete functionFactory;
258 return;
259 }
260
261 AIDA::IPlotter* plotter = plotterFactory->create();
262 if(!plotter) {
263 Lib::Out out(fSession.printer());
264 out << "Can't create a plotter factory." << Lib::endl;
265 delete functionFactory;
266 delete plotterFactory;
267 return;
268 }
269
270 int number = obj->fTable->entries();
271 if(!superpose) {
272 if(number>0) plotter->createRegions(1,number,0);
273 else plotter->createRegions(1,1,0);
274 }
275
276 for (int count=0;count<number;count++){
277 G4PhysicsVector* physicsVector = (*(obj->fTable))[count];
278 size_t n = physicsVector->GetVectorLength();
279 if(n>=2) {
280 std::string name;
281 Lib::smanip::printf(name,1024,"%s.%s.%s.%d",
282 obj->fTableName.c_str(),
283 obj->fProcessName.c_str(),
284 obj->fParticleName.c_str(),
285 count);
286 std::vector<double> params;
287 params.push_back(1); // Dimension of the grid (here one).
288 params.push_back(n); // Number of entries for first grid axis.
289 size_t i;
290 for (i=0;i<n;i++) params.push_back(i); // Xs.
291 for (i=0;i<n;i++) params.push_back((*physicsVector)[i]); //Values.
292 AIDA::IFunction* function =
293 functionFactory->createFunctionByName(name,"Grid1D");
294 if(!function) {
295 Lib::Out out(fSession.printer());
296 out << "Can't create \"" << name << "\" function." << Lib::endl;
297 break;
298 }
299 function->setParameters(params);
300
301 AIDA::IPlotterRegion& region = plotter->currentRegion();
302
303 region.plot(*function);
304 region.setParameter("plotter.xAxisAutomated","FALSE");
305 std::string value1;
306 Lib::smanip::printf(value1,32,"%g",0.);
307 region.setParameter("plotter.xAxisMinimum",value1);
308 std::string value2;
309 Lib::smanip::printf(value2,32,"%g",(double)(n-1));
310 region.setParameter("plotter.xAxisMaximum",value2);
311 std::string value3;
312 Lib::smanip::printf(value3,32,"%d",4 * n);
313 region.setParameter("plotter.xNumberOfPoints",value3);
314 /*
315 fSession.out().println("debug : min : %s, max : %s, steps : %s\n.",
316 value1.c_str(),value2.c_str(),value3.c_str());
317 */
318
319 if(!superpose) plotter->next();
320
321 }
322 }
323
324 delete plotter;
325 delete functionFactory;
326 delete plotterFactory;
327
328}
329//////////////////////////////////////////////////////////////////////////////
330void G4Lab::PhysicsTableAccessor::endVisualize(
331)
332//////////////////////////////////////////////////////////////////////////////
333//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
334{
335}
336//////////////////////////////////////////////////////////////////////////////
337bool G4Lab::PhysicsTableAccessor::findProcessTables(
338 G4ParticleDefinition* aParticle
339,G4VProcess* aProcess
340,std::vector<std::string>& aTables
341,std::vector<G4PhysicsTable*>& aPhysicsTables
342)
343//////////////////////////////////////////////////////////////////////////////
344//!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!//
345{
346 aTables.clear();
347 aPhysicsTables.clear();
348
349 std::string tmpName;
350 if(!Lib::fmanip::tmpname(tmpName,".","G4Lab_tmp","")) {
351 Lib::Out out(fSession.printer());
352 out << "Can't build a temporary directory name." << Lib::endl;
353 return false;
354 }
355
356 if(!Lib::dirmanip::create(tmpName)) {
357 Lib::Out out(fSession.printer());
358 out << "Can't create \"" << tmpName << "\" directory." << Lib::endl;
359 return false;
360 }
361
362#ifdef WIN32 //binary mode out on Windows.
363 G4bool ascii = true;
364#else
365 G4bool ascii = false;
366#endif
367 aProcess->StorePhysicsTable(aParticle,tmpName,ascii);
368
369 // Analyse the content of the directory :
370 std::vector<std::string> files;
371 if(!Lib::dirmanip::entries(tmpName,files,false)) {
372 Lib::Out out(fSession.printer());
373 out << "Can't get files in \"" << tmpName << "\" directory." << Lib::endl;
374 } else {
375 // Analyse files :
376 Lib::smanip::remove(files,".");
377 Lib::smanip::remove(files,"..");
378 for(unsigned int index=0;index<files.size();index++) {
379 std::vector<std::string> words;
380 Lib::smanip::words(files[index],".",words);
381 //Lib::Out out(fSession.printer());
382 //out << "debug : Analyse \"" << files[index] << "\" ." << Lib::endl;
383 if(words.size()==4) {
384 std::string table = words[0];
385 std::string process = words[1];
386 std::string particle = words[2];
387 //Lib::smanip::replace(aTables[i],".dat","");
388 std::string filename =
389 aProcess->GetPhysicsTableFileName(aParticle,
390 tmpName,
391 table,
392 ascii);
393 G4PhysicsTable* physicsTable =
394 new G4PhysicsTable(G4Element::GetNumberOfElements());
395 if(physicsTable) {
396 G4bool status = physicsTable->RetrievePhysicsTable(filename,ascii);
397 if(!status){
398 Lib::Out out(fSession.printer());
399 out << "Unable to read \"" << filename << "\" file." << Lib::endl;
400 delete physicsTable;
401 physicsTable = 0;
402 }
403 }
404 if(physicsTable) {
405 aTables.push_back(table+"."+process+"."+particle);
406 aPhysicsTables.push_back(physicsTable);
407 }
408 }
409 }
410 }
411
412 // Remove the temporary directory :
413 if(!Lib::dirmanip::remove(tmpName)) {
414 Lib::Out out(fSession.printer());
415 out << "Can't remove \"" << tmpName << "\" directory." << Lib::endl;
416 aTables.clear();
417 aPhysicsTables.clear();
418 return false;
419 }
420
421 return true;
422}
Note: See TracBrowser for help on using the repository browser.