Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TFoamCell.cxx
Go to the documentation of this file.
1// @(#)root/foam:$Id$
2// Author: S. Jadach <mailto:Stanislaw.jadach@ifj.edu.pl>, P.Sawicki <mailto:Pawel.Sawicki@ifj.edu.pl>
3
4/** \class TFoamCell
5
6Used by TFoam
7
8Objects of this class are hyper-rectangular cells organized in the binary tree.
9Special algorithm for encoding relative positioning of the cells
10allow to save total memory allocation needed for the system of cells.
11*/
12
13
14#include <iostream>
15#include "TFoamCell.h"
16#include "TFoamVect.h"
17
18
20
21////////////////////////////////////////////////////////////////////////////////
22/// Default constructor for streamer
23
25{
26 fParent = 0;
27 fDaught0 = 0;
28 fDaught1 = 0;
29}
30
31////////////////////////////////////////////////////////////////////////////////
32/// User constructor allocating single empty Cell
33
35{
36 if ( kDim >0) {
37 //---------=========----------
38 fDim = kDim;
39 fSerial = 0;
40 fStatus = 1;
41 fParent = 0;
42 fDaught0 = 0;
43 fDaught1 = 0;
44 fXdiv = 0.0;
45 fBest = 0;
46 fVolume = 0.0;
47 fIntegral = 0.0;
48 fDrive = 0.0;
49 fPrimary = 0.0;
50 } else
51 Error("TFoamCell","Dimension has to be >0 \n ");
52}
53
54////////////////////////////////////////////////////////////////////////////////
55/// Copy constructor (not tested!)
56
58{
59 Error("TFoamCell", "+++++ NEVER USE Copy constructor for TFoamCell \n");
60 fStatus = From.fStatus;
61 fParent = From.fParent;
62 fDaught0 = From.fDaught0;
63 fDaught1 = From.fDaught1;
64 fXdiv = From.fXdiv;
65 fBest = From.fBest;
66 fVolume = From.fVolume;
67 fIntegral = From.fIntegral;
68 fDrive = From.fDrive;
69 fPrimary = From.fPrimary;
70}
71
72////////////////////////////////////////////////////////////////////////////////
73/// Destructor
74
76{
77}
78
79////////////////////////////////////////////////////////////////////////////////
80/// Substitution operator = (never used)
81
83{
84 Info("TFoamCell", "operator=\n ");
85 if (&From == this) return *this;
86 fStatus = From.fStatus;
87 fParent = From.fParent;
88 fDaught0 = From.fDaught0;
89 fDaught1 = From.fDaught1;
90 fXdiv = From.fXdiv;
91 fBest = From.fBest;
92 fVolume = From.fVolume;
93 fIntegral = From.fIntegral;
94 fDrive = From.fDrive;
95 fPrimary = From.fPrimary;
96 return *this;
97}
98
99
100////////////////////////////////////////////////////////////////////////////////
101/// Fills in certain data into newly allocated cell
102
103void TFoamCell::Fill(Int_t Status, TFoamCell *Parent, TFoamCell *Daugh1, TFoamCell *Daugh2)
104{
105 fStatus = Status;
106 fParent = Parent;
107 fDaught0 = Daugh1;
108 fDaught1 = Daugh2;
109}
110
111////////////////////////////////////////////////////////////////////////////////
112/// Provides size and position of the cell
113/// These parameter are calculated by analyzing information in all parents
114/// cells up to the root cell. It takes time but saves memory.
115
116void TFoamCell::GetHcub( TFoamVect &cellPosi, TFoamVect &cellSize) const
117{
118 if(fDim<1) return;
119 const TFoamCell *pCell,*dCell;
120 cellPosi = 0.0; cellSize=1.0; // load all components
121 dCell = this;
122 while(dCell != 0) {
123 pCell = dCell->GetPare();
124 if( pCell== 0) break;
125 Int_t kDiv = pCell->fBest;
126 Double_t xDivi = pCell->fXdiv;
127 if(dCell == pCell->GetDau0() ) {
128 cellSize[kDiv] *=xDivi;
129 cellPosi[kDiv] *=xDivi;
130 } else if( dCell == pCell->GetDau1() ) {
131 cellSize[kDiv] *=(1.0-xDivi);
132 cellPosi[kDiv] =cellPosi[kDiv]*(1.0-xDivi)+xDivi;
133 } else {
134 Error("GetHcub ","Something wrong with linked tree \n");
135 }
136 dCell=pCell;
137 }//while
138}
139
140////////////////////////////////////////////////////////////////////////////////
141/// Provides size of the cell
142/// Size parameters are calculated by analyzing information in all parents
143/// cells up to the root cell. It takes time but saves memory.
144
145void TFoamCell::GetHSize( TFoamVect &cellSize) const
146{
147 if(fDim<1) return;
148 const TFoamCell *pCell,*dCell;
149 cellSize=1.0; // load all components
150 dCell = this;
151 while(dCell != 0) {
152 pCell = dCell->GetPare();
153 if( pCell== 0) break;
154 Int_t kDiv = pCell->fBest;
155 Double_t xDivi = pCell->fXdiv;
156 if(dCell == pCell->GetDau0() ) {
157 cellSize[kDiv]=cellSize[kDiv]*xDivi;
158 } else if(dCell == pCell->GetDau1() ) {
159 cellSize[kDiv]=cellSize[kDiv]*(1.0-xDivi);
160 } else {
161 Error("GetHSize ","Something wrong with linked tree \n");
162 }
163 dCell=pCell;
164 }//while
165}
166
167////////////////////////////////////////////////////////////////////////////////
168/// Calculates volume of the cell using size params which are calculated
169
171{
172 Int_t k;
173 Double_t volu=1.0;
174 if(fDim>0) { // h-cubical subspace
175 TFoamVect cellSize(fDim);
176 GetHSize(cellSize);
177 for(k=0; k<fDim; k++) volu *= cellSize[k];
178 }
179 fVolume =volu;
180}
181
182////////////////////////////////////////////////////////////////////////////////
183/// Printout of the cell geometry parameters for the debug purpose
184
185void TFoamCell::Print(Option_t *option) const
186{
187 if(!option) Error("Print", "No option set\n");
188
189 std::cout << " Status= "<< fStatus <<",";
190 std::cout << " Volume= "<< fVolume <<",";
191 std::cout << " TrueInteg= " << fIntegral <<",";
192 std::cout << " DriveInteg= "<< fDrive <<",";
193 std::cout << " PrimInteg= " << fPrimary <<",";
194 std::cout<< std::endl;
195 std::cout << " Xdiv= "<<fXdiv<<",";
196 std::cout << " Best= "<<fBest<<",";
197 std::cout << " Parent= {"<< (GetPare() ? GetPare()->GetSerial() : -1) <<"} "; // extra DEBUG
198 std::cout << " Daught0= {"<< (GetDau0() ? GetDau0()->GetSerial() : -1 )<<"} "; // extra DEBUG
199 std::cout << " Daught1= {"<< (GetDau1() ? GetDau1()->GetSerial() : -1 )<<"} "; // extra DEBUG
200 std::cout<< std::endl;
201 //
202 //
203 if(fDim>0 ) {
204 TFoamVect cellPosi(fDim); TFoamVect cellSize(fDim);
205 GetHcub(cellPosi,cellSize);
206 std::cout <<" Posi= "; cellPosi.Print("1"); std::cout<<","<< std::endl;
207 std::cout <<" Size= "; cellSize.Print("1"); std::cout<<","<< std::endl;
208 }
209}
double Double_t
Definition RtypesCore.h:59
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:364
const Int_t kDiv
Definition TPCON.h:30
Used by TFoam.
Definition TFoamCell.h:12
Double_t fXdiv
Factor for division.
Definition TFoamCell.h:26
virtual ~TFoamCell()
Destructor.
Definition TFoamCell.cxx:75
TFoamCell * GetDau1() const
Definition TFoamCell.h:62
void CalcVolume()
Calculates volume of the cell using size params which are calculated.
void Print(Option_t *option) const
Printout of the cell geometry parameters for the debug purpose.
TRef fDaught1
Pointer to daughter 2.
Definition TFoamCell.h:23
Int_t GetSerial() const
Definition TFoamCell.h:66
Double_t fIntegral
Integral over cell (estimate from exploration)
Definition TFoamCell.h:30
Int_t fStatus
Status (active, inactive)
Definition TFoamCell.h:20
TFoamCell * GetDau0() const
Definition TFoamCell.h:61
Double_t fVolume
Cartesian Volume of cell.
Definition TFoamCell.h:29
TRef fParent
Pointer to parent cell.
Definition TFoamCell.h:21
TRef fDaught0
Pointer to daughter 1.
Definition TFoamCell.h:22
void GetHcub(TFoamVect &, TFoamVect &) const
Provides size and position of the cell These parameter are calculated by analyzing information in all...
TFoamCell * GetPare() const
Definition TFoamCell.h:60
Double_t fDrive
Driver integral, only for cell build-up.
Definition TFoamCell.h:31
Short_t fDim
Dimension of the vector space.
Definition TFoamCell.h:15
Double_t fPrimary
Primary integral, only for MC generation.
Definition TFoamCell.h:32
TFoamCell()
Default constructor for streamer.
Definition TFoamCell.cxx:24
TFoamCell & operator=(const TFoamCell &)
Substitution operator = (never used)
Definition TFoamCell.cxx:82
void Fill(Int_t, TFoamCell *, TFoamCell *, TFoamCell *)
Fills in certain data into newly allocated cell.
Int_t fBest
Best Edge for division.
Definition TFoamCell.h:27
Int_t fSerial
Serial number.
Definition TFoamCell.h:19
void GetHSize(TFoamVect &) const
Provides size of the cell Size parameters are calculated by analyzing information in all parents cell...
Auxiliary class TFoamVect of n-dimensional vector, with dynamic allocation used for the cartesian geo...
Definition TFoamVect.h:10
void Print(Option_t *option) const
Printout of all vector components on "std::cout".
Mother of all ROOT objects.
Definition TObject.h:41
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:963
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:937