Jpp test-rotations-old
the software that should make you happy
Loading...
Searching...
No Matches
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
JMARKOV::JMarkovEnsembleIntegrator1D Class Reference

Implementation of JMarkovEnsembleIntegrator interface with 1D histograms. More...

#include <JMarkovIntegrator.hh>

Inheritance diagram for JMARKOV::JMarkovEnsembleIntegrator1D:
JMARKOV::JMarkovEnsembleIntegrator JMARKOV::JMarkovIntegrator

Public Member Functions

 JMarkovEnsembleIntegrator1D (const vector< JPhotonPath > &_ensemble, JTargetModel *_trg, int _nbinsx, int _nbinsy, int _nbinsz)
 
 ~JMarkovEnsembleIntegrator1D ()
 
void free_scattering_vertex_gens ()
 
void writeHistograms ()
 write the histograms filled from the ensemble
 
vector< double > integrate (int N, int nscat, JSourceModel *src, JScatteringModel *sm, JTargetModel *trg, double lambda_abs)
 Integrate with N samples.
 
vector< double > dummy_integrate (int N, int nscat, JSourceModel *src, JScatteringModel *sm, JTargetModel *trg, double lambda_abs)
 Integrate a test function with N samples.
 
vector< JPhotonPath > get_diagnostic_ensemble (int N, int nscat, JSourceModel *src, JScatteringModel *sm, JTargetModel *trg, double lambda_abs)
 Return photon paths generated with the generatePath method.
 

Protected Member Functions

void init_scattering_vertex_generators (int nscat)
 
JPosition3D generateScatteringVertexPosition (int nscat, int nv, double &winv)
 
JPosition3D generatePosition (int nscat, int nv, double &winv)
 generate a position for vertex nv for the given number of scatterings, winv will be set to 1/weight
 
int getIndex (int nv, int nscat)
 return the internal index for a given number of scatterings and vertex number
 
void initgenerators (int nscat)
 initialize the position generators for a given number of scatterings (i.e. create histograms and fill them from the ensemble)
 
void init_target_generator (int nscat)
 initialize the generator for a target vertex for a given number of scatterings
 
void free_target_gens ()
 
virtual JPhotonPath generatePath (int nscat, double &winv)
 Generate a random photon path with a given number of scatterings.
 

Protected Attributes

vector< JHistGenerator * > gens
 
vector< JSphereGenerator * > target_gens
 
vector< JPhotonPath > ensemble
 
int nbinsx
 
int nbinsy
 
int nbinsz
 
JTargetModel * trg
 
JPosition3D source_position
 

Detailed Description

Implementation of JMarkovEnsembleIntegrator interface with 1D histograms.

Definition at line 313 of file JMarkovIntegrator.hh.

Constructor & Destructor Documentation

◆ JMarkovEnsembleIntegrator1D()

JMARKOV::JMarkovEnsembleIntegrator1D::JMarkovEnsembleIntegrator1D ( const vector< JPhotonPath > & _ensemble,
JTargetModel * _trg,
int _nbinsx,
int _nbinsy,
int _nbinsz )
inline

Definition at line 317 of file JMarkovIntegrator.hh.

317: JMarkovEnsembleIntegrator(_ensemble,_trg,_nbinsx,_nbinsy,_nbinsz) {}
JMARKOV::JMarkovEnsembleIntegrator::JMarkovEnsembleIntegrator
JMarkovEnsembleIntegrator(const vector< JPhotonPath > &_ensemble, JTargetModel *_trg, int _nbinsx, int _nbinsy, int _nbinsz)
Definition JMarkovIntegrator.hh:202

◆ ~JMarkovEnsembleIntegrator1D()

JMARKOV::JMarkovEnsembleIntegrator1D::~JMarkovEnsembleIntegrator1D ( )
inline

Definition at line 319 of file JMarkovIntegrator.hh.

319 {
320 free_target_gens() ;
321 free_scattering_vertex_gens() ;
322 }

Member Function Documentation

◆ free_scattering_vertex_gens()

void JMARKOV::JMarkovEnsembleIntegrator1D::free_scattering_vertex_gens ( )
inline

Definition at line 324 of file JMarkovIntegrator.hh.

324 {
325 for( vector<JHistGenerator*>::iterator it=gens.begin() ; it!=gens.end() ; ++it ) {
326 if( *it != NULL ) delete *it ;
327 }
328 }

◆ writeHistograms()

void JMARKOV::JMarkovEnsembleIntegrator1D::writeHistograms ( )
inlinevirtual

write the histograms filled from the ensemble

Implements JMARKOV::JMarkovEnsembleIntegrator.

Definition at line 330 of file JMarkovIntegrator.hh.

330 {
331 for( vector<JHistGenerator*>::iterator it=gens.begin() ; it!=gens.end() ; ++it ) {
332 if( *it != NULL ) {
333 (*it)->hx->Write() ;
334 (*it)->hy->Write() ;
335 (*it)->hz->Write() ;
336 }
337 }
338 for( vector<JSphereGenerator*>::iterator it=target_gens.begin(); it!=target_gens.end(); ++it ) {
339 if( *it != NULL ) {
340 if( (*it)->h != NULL ) (*it)->h->Write() ;
341 }
342 }
343 }
JMARKOV::JMarkovEnsembleIntegrator::target_gens
vector< JSphereGenerator * > target_gens
Definition JMarkovIntegrator.hh:301

◆ init_scattering_vertex_generators()

void JMARKOV::JMarkovEnsembleIntegrator1D::init_scattering_vertex_generators ( int nscat)
inlineprotectedvirtual

Implements JMARKOV::JMarkovEnsembleIntegrator.

Definition at line 347 of file JMarkovIntegrator.hh.

347 {
348 int ifirst = getIndex(1,nscat) ;
349 int ilast = getIndex(1,nscat+1) ;
350 if( (int)gens.size() < ilast ) gens.resize( ilast, NULL ) ;
351 for( int nv=1 ; nv<=nscat ; ++nv ) {
352 int index = ifirst + nv - 1 ;
353 // determine minimal and maximal coordinate values for this vertex
354 double xmin = 1.0/0.0 ;
355 double xmax = -1.0/0.0 ;
356 double ymin = 1.0/0.0 ;
357 double ymax = -1.0/0.0 ;
358 double zmin = 1.0/0.0 ;
359 double zmax = -1.0/0.0 ;
360 for( vector<JPhotonPath>::iterator it=ensemble.begin() ; it!=ensemble.end() ; ++it ) {
361 if( it->n-2 != nscat ) continue ;
362 JPosition3D vtx( (*it)[nv] ) ;
363 if( vtx.getX()<xmin ) xmin = vtx.getX() ;
364 if( vtx.getX()>xmax ) xmax = vtx.getX() ;
365 if( vtx.getY()<ymin ) ymin = vtx.getY() ;
366 if( vtx.getY()>ymax ) ymax = vtx.getY() ;
367 if( vtx.getZ()<zmin ) zmin = vtx.getZ() ;
368 if( vtx.getZ()>zmax ) zmax = vtx.getZ() ;
369 }
370 // allocate histograms
371 char hname[200] ;
372 sprintf( hname, "hx_nscat%i_vtx%i", nscat, nv ) ;
373 TH1F hx(hname,"Vertex position;X",nbinsx,xmin,xmax) ;
374 sprintf( hname, "hy_nscat%i_vtx%i", nscat, nv ) ;
375 TH1F hy(hname,"Vertex position;Y",nbinsy,ymin,ymax) ;
376 sprintf( hname, "hz_nscat%i_vtx%i", nscat, nv ) ;
377 TH1F hz(hname,"Vertex position;Z",nbinsz,zmin,zmax) ;
378 // fill histograms from ensemble
379 for( vector<JPhotonPath>::iterator it=ensemble.begin() ; it!=ensemble.end() ; ++it ) {
380 if( it->n-2 != nscat ) continue ;
381 JPosition3D vtx( (*it)[nv] ) ;
382 hx.Fill( vtx.getX() ) ;
383 hy.Fill( vtx.getY() ) ;
384 hz.Fill( vtx.getZ() ) ;
385 }
386 // extra: make sure every field is filled at least once
387 for( Int_t bin=1 ; bin<=hx.GetNbinsX() ; ++bin )
388 hx.AddBinContent(bin,1) ;
389 for( Int_t bin=1 ; bin<=hy.GetNbinsX() ; ++bin )
390 hy.AddBinContent(bin,1) ;
391 for( Int_t bin=1 ; bin<=hz.GetNbinsX() ; ++bin )
392 hz.AddBinContent(bin,1) ;
393
394 gens[index] = new JHistGenerator( &hx, &hy, &hz ) ;
395 }
396 }
JMARKOV::JMarkovEnsembleIntegrator::getIndex
int getIndex(int nv, int nscat)
return the internal index for a given number of scatterings and vertex number
Definition JMarkovIntegrator.hh:251
GAUSS_LEGENDRE::xmax
const double xmax
Definition JQuadrature.cc:24
GAUSS_LEGENDRE::xmin
const double xmin
Definition JQuadrature.cc:23

◆ generateScatteringVertexPosition()

JPosition3D JMARKOV::JMarkovEnsembleIntegrator1D::generateScatteringVertexPosition ( int nscat,
int nv,
double & winv )
inlineprotectedvirtual

Implements JMARKOV::JMarkovEnsembleIntegrator.

Definition at line 399 of file JMarkovIntegrator.hh.

399 {
400 int index = getIndex(nv,nscat) ;
401 if( index+1>(int)gens.size() || gens[index] == NULL ) initgenerators(nscat) ;
402
403 JPosition3D pos = gens[index]->getPosition() ;
404 winv = 1.0/gens[index]->getWeight(pos) ;
405 return pos ;
406 }
JMARKOV::JMarkovEnsembleIntegrator::initgenerators
void initgenerators(int nscat)
initialize the position generators for a given number of scatterings (i.e. create histograms and fill...
Definition JMarkovIntegrator.hh:263

◆ generatePosition()

JPosition3D JMARKOV::JMarkovEnsembleIntegrator::generatePosition ( int nscat,
int nv,
double & winv )
inlineprotectedvirtualinherited

generate a position for vertex nv for the given number of scatterings, winv will be set to 1/weight

Implements JMARKOV::JMarkovIntegrator.

Definition at line 228 of file JMarkovIntegrator.hh.

228 {
229 // source vertex
230 if( nv==0 ) {
231 winv = 1.0 ;
232 return source_position ;
233 }
234 // target vertex
235 if( nv == nscat+1 ) {
236 if( nscat+1>(int)target_gens.size() || target_gens[nscat] == NULL ) {
237 initgenerators(nscat) ;
238 }
239 JPosition3D pos = target_gens[nscat]->getPosition() ;
240 winv = 1.0/target_gens[nscat]->getWeight(pos) ;
241 return pos ;
242 } else {
243 // scattering vertex
244 return generateScatteringVertexPosition(nscat,nv,winv) ;
245 }
246 }
JMARKOV::JMarkovEnsembleIntegrator::generateScatteringVertexPosition
virtual JPosition3D generateScatteringVertexPosition(int nscat, int nv, double &winv)=0

◆ getIndex()

int JMARKOV::JMarkovEnsembleIntegrator::getIndex ( int nv,
int nscat )
inlineprotectedinherited

return the internal index for a given number of scatterings and vertex number

Definition at line 251 of file JMarkovIntegrator.hh.

251 {
252 // indices are distributed as follows
253 // [0]: nscat=1, nv=1
254 // [1]: nscat=2, nv=1
255 // [2]: nscat=2, nv=2
256 // [3]: nscat=3, nv=1
257 // [4]: nscat=3, nv=2
258 // etc.
259 return (nscat*(nscat-1))/2 + nv - 1 ;
260 }

◆ initgenerators()

void JMARKOV::JMarkovEnsembleIntegrator::initgenerators ( int nscat)
inlineprotectedinherited

initialize the position generators for a given number of scatterings (i.e. create histograms and fill them from the ensemble)

Definition at line 263 of file JMarkovIntegrator.hh.

263 {
264 // target vertex
265 init_target_generator(nscat) ;
266
267 // scattering vertices
268 init_scattering_vertex_generators(nscat) ;
269 }
JMARKOV::JMarkovEnsembleIntegrator::init_scattering_vertex_generators
virtual void init_scattering_vertex_generators(int nscat)=0
JMARKOV::JMarkovEnsembleIntegrator::init_target_generator
void init_target_generator(int nscat)
initialize the generator for a target vertex for a given number of scatterings
Definition JMarkovIntegrator.hh:272

◆ init_target_generator()

void JMARKOV::JMarkovEnsembleIntegrator::init_target_generator ( int nscat)
inlineprotectedinherited

initialize the generator for a target vertex for a given number of scatterings

Definition at line 272 of file JMarkovIntegrator.hh.

272 {
273 if( (int)target_gens.size() <= nscat ) target_gens.resize( nscat+1 ) ;
274
275 if( trg->getRadius()>0 ) {
276 // for a finite target, we fill a "sphere generator" from the ensemble
277 char hname[200] ;
278 sprintf( hname, "htarget_nscat%i", nscat ) ;
279 TH2D* ht = new TH2D(hname,";cos(#theta);#phi",100,-1,1,100,-M_PI,M_PI) ;
280 for( vector<JPhotonPath>::iterator it=ensemble.begin() ; it!=ensemble.end() ; ++it ) {
281 if( it->n-2 != nscat ) continue ;
282 JDirection3D dir( it->back() - trg->getPosition() ) ;
283 ht->Fill( dir.getDZ(), dir.getPhi() ) ;
284 }
285 target_gens[nscat] = new JSphereGenerator( trg->getPosition(), trg->getRadius(), ht ) ;
286 delete ht ;
287 } else {
288 // for an infinitesimal target, we create a "generator" that only generates the single point
289 target_gens[nscat] = new JSphereGenerator( trg->getPosition() ) ;
290 }
291 }
JMARKOV::JTargetModel::getPosition
const JPosition3D & getPosition() const
Definition JScatteringModel.hh:304

◆ free_target_gens()

void JMARKOV::JMarkovEnsembleIntegrator::free_target_gens ( )
inlineprotectedinherited

Definition at line 295 of file JMarkovIntegrator.hh.

295 {
296 for( vector<JSphereGenerator*>::iterator it=target_gens.begin() ; it!=target_gens.end() ; ++it ) {
297 if( *it != NULL ) delete *it ;
298 }
299 }

◆ integrate()

vector< double > JMARKOV::JMarkovIntegrator::integrate ( int N,
int nscat,
JSourceModel * src,
JScatteringModel * sm,
JTargetModel * trg,
double lambda_abs )
inherited

Integrate with N samples.

Returns a vector with the contribution to the integral of each sample. The mean of those values is the estimate of the result of the integral, while the distribution itself can be used to estimate the stability of the result. In this distribution, you want to avoid

  • long tails (because they make the result unstable)
  • small contributions (because it means that parts of the parameter space are being oversampled, so it is less efficient) This can be achieved by tuning the sample distribution to the problem at hand.

Definition at line 130 of file JMarkovIntegrator.hh.

130 {
131 vector<double> contributions(N,-1) ;
132
133 for( int i=0 ; i<N ; ++i ) {
134 double winv ;
135 JPhotonPath p = generatePath(nscat,winv) ;
136
137 double rho = getPhotonPathProbabilityDensity(p,src,sm,trg,lambda_abs) ;
138 contributions[i] = rho * winv ;
139 }
140 return contributions ;
141 }
JMARKOV::JMarkovIntegrator::generatePath
virtual JPhotonPath generatePath(int nscat, double &winv)
Generate a random photon path with a given number of scatterings.
Definition JMarkovIntegrator.hh:94
JMARKOV::getPhotonPathProbabilityDensity
double getPhotonPathProbabilityDensity(JPhotonPath &p, JSourceModel *src, JScatteringModel *sm, JTargetModel *trg, double lambda_abs)
Return the probability density for a photon path with the given ingredients.
Definition JScatteringModel.hh:582

◆ dummy_integrate()

vector< double > JMARKOV::JMarkovIntegrator::dummy_integrate ( int N,
int nscat,
JSourceModel * src,
JScatteringModel * sm,
JTargetModel * trg,
double lambda_abs )
inherited

Integrate a test function with N samples.

This can be used as a sanity check for derived classes of JMarkovIntegrator.

The integral should yield 1 when the complete relevant part of the volume is taken into account. If it does not, it may be a sign that the implementation is not correct.

Returns a vector with the contribution to the integral of each sample.

Definition at line 143 of file JMarkovIntegrator.hh.

143 {
144 vector<double> contributions(N,-1) ;
145
146 for( int i=0 ; i<N ; ++i ) {
147 const double r = 10 ;
148 double winv ;
149 JPhotonPath p = generatePath(nscat,winv) ;
150 double rho = 1.0/(4.0/3.0*M_PI*r*r*r) ;
151 if( p[1].getLength()>r ) rho = 0 ;
152 contributions[i] = rho * winv ;
153 }
154 return contributions ;
155 }

◆ get_diagnostic_ensemble()

vector< JPhotonPath > JMARKOV::JMarkovIntegrator::get_diagnostic_ensemble ( int N,
int nscat,
JSourceModel * src,
JScatteringModel * sm,
JTargetModel * trg,
double lambda_abs )
inherited

Return photon paths generated with the generatePath method.

This can be used to identify the parts of parameter space that are over- or undersampled in a given problem so that the integrator may be optimized to handle those better.

Definition at line 157 of file JMarkovIntegrator.hh.

157 {
158 vector<JPhotonPath> paths ;
159
160 for( int i=0 ; i<N ; ++i ) {
161 double winv ;
162 JPhotonPath p = generatePath(nscat,winv) ;
163 paths.push_back(p) ;
164 }
165 return paths ;
166 }

◆ generatePath()

virtual JPhotonPath JMARKOV::JMarkovIntegrator::generatePath ( int nscat,
double & winv )
inlineprotectedvirtualinherited

Generate a random photon path with a given number of scatterings.

winv must be set to the inverted probability density to generate this particular path.

Definition at line 94 of file JMarkovIntegrator.hh.

94 {
95 // this default implementation assumes that the vertex positions
96 // are completely uncorrelated
97 JPhotonPath p(nscat) ;
98 double _winv = 1 ;
99 for( int nv=0 ; nv<nscat+2 ; ++nv ) {
100 double part_winv ;
101 p[nv] = generatePosition( nscat, nv, part_winv ) ;
102 _winv *= part_winv ;
103 }
104 winv = _winv ;
105 return p ;
106 }
JMARKOV::JMarkovIntegrator::generatePosition
virtual JPosition3D generatePosition(int nscat, int nv, double &winv)=0
Generate a random position for vertex nv.

Member Data Documentation

◆ gens

vector<JHistGenerator*> JMARKOV::JMarkovEnsembleIntegrator1D::gens
protected

Definition at line 408 of file JMarkovIntegrator.hh.

◆ target_gens

vector<JSphereGenerator*> JMARKOV::JMarkovEnsembleIntegrator::target_gens
protectedinherited

Definition at line 301 of file JMarkovIntegrator.hh.

◆ ensemble

vector<JPhotonPath> JMARKOV::JMarkovEnsembleIntegrator::ensemble
protectedinherited

Definition at line 302 of file JMarkovIntegrator.hh.

◆ nbinsx

int JMARKOV::JMarkovEnsembleIntegrator::nbinsx
protectedinherited

Definition at line 303 of file JMarkovIntegrator.hh.

◆ nbinsy

int JMARKOV::JMarkovEnsembleIntegrator::nbinsy
protectedinherited

Definition at line 304 of file JMarkovIntegrator.hh.

◆ nbinsz

int JMARKOV::JMarkovEnsembleIntegrator::nbinsz
protectedinherited

Definition at line 305 of file JMarkovIntegrator.hh.

◆ trg

JTargetModel* JMARKOV::JMarkovEnsembleIntegrator::trg
protectedinherited

Definition at line 306 of file JMarkovIntegrator.hh.

◆ source_position

JPosition3D JMARKOV::JMarkovEnsembleIntegrator::source_position
protectedinherited

Definition at line 307 of file JMarkovIntegrator.hh.

The documentation for this class was generated from the following file:
Generated by doxygen 1.12.0