Abstract base class for implementations of the JMarkovIntegrator interface, where the sample distribution is based on histograms filled from an ensemble of representative paths. More...

#include <JMarkovIntegrator.hh>

Inheritance diagram for JMARKOV::JMarkovEnsembleIntegrator:

Public Member Functions
	JMarkovEnsembleIntegrator (const vector< JPhotonPath > &_ensemble, JTargetModel *_trg, int _nbinsx, int _nbinsy, int _nbinsz)

virtual	~JMarkovEnsembleIntegrator ()

virtual void	writeHistograms ()=0
	write the histograms filled from the ensemble More...

vector< double >	integrate (int N, int nscat, JSourceModel src, JScatteringModel sm, JTargetModel *trg, double lambda_abs)
	Integrate with N samples. More...

vector< double >	dummy_integrate (int N, int nscat, JSourceModel src, JScatteringModel sm, JTargetModel *trg, double lambda_abs)
	Integrate a test function with N samples. More...

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. More...

Protected Member Functions
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 More...

virtual JPosition3D	generateScatteringVertexPosition (int nscat, int nv, double &winv)=0

int	getIndex (int nv, int nscat)
	return the internal index for a given number of scatterings and vertex number More...

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

void	init_target_generator (int nscat)
	initialize the generator for a target vertex for a given number of scatterings More...

virtual void	init_scattering_vertex_generators (int nscat)=0

void	free_target_gens ()

virtual JPhotonPath	generatePath (int nscat, double &winv)
	Generate a random photon path with a given number of scatterings. More...

Protected Attributes
vector< JSphereGenerator * >	target_gens

vector< JPhotonPath >	ensemble

int	nbinsx

int	nbinsy

int	nbinsz

JTargetModel *	trg

JPosition3D	source_position

Detailed Description

Abstract base class for implementations of the JMarkovIntegrator interface, where the sample distribution is based on histograms filled from an ensemble of representative paths.

Definition at line 198 of file JMarkovIntegrator.hh.

Constructor & Destructor Documentation

◆ JMarkovEnsembleIntegrator()

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

inline

Definition at line 202 of file JMarkovIntegrator.hh.

                                                                                                                                  : ensemble(_ensemble), nbinsx(_nbinsx), nbinsy(_nbinsy), nbinsz(_nbinsz), trg(_trg) {
       // check that the ensemble size is at least one photon path
       if( ensemble.size()==0 ) {
         cerr << "FATAL ERROR in JMarkov3DensembleIntegrator constructor. Ensemble size = 0." << endl ;
         exit(1) ;
       }
  
       // check that the position of the first vertex is the same throughout the ensemble
       for( vector<JPhotonPath>::iterator it=ensemble.begin(); it!=ensemble.end(); ++it ) {
         if( (*it)[0].getDistance(ensemble.front()[0]) > 0 ) {
           cerr << "Fatal error in JMarkovEnsembleIntegrator, first vertex positions are not all the same." << endl ;
           exit(1) ;
         }
       }
       // set source position based on ensemble
       source_position = ensemble.front()[0] ;
     }

◆ ~JMarkovEnsembleIntegrator()

virtual JMARKOV::JMarkovEnsembleIntegrator::~JMarkovEnsembleIntegrator ( )

inlinevirtual

Definition at line 220 of file JMarkovIntegrator.hh.

220 {} ;

Member Function Documentation

◆ writeHistograms()

virtual void JMARKOV::JMarkovEnsembleIntegrator::writeHistograms ( )

pure virtual

write the histograms filled from the ensemble

Implemented in JMARKOV::JMarkovEnsembleIntegrator3D, and JMARKOV::JMarkovEnsembleIntegrator1D.

◆ generatePosition()

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

inlineprotectedvirtual

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.

                                                                     {
       // source vertex
       if( nv==0 ) {
         winv = 1.0 ;
         return source_position ;
       }
       // target vertex
       if( nv == nscat+1 ) {
         if( nscat+1>(int)target_gens.size() || target_gens[nscat] == NULL ) { 
           initgenerators(nscat) ; 
         }
         JPosition3D pos = target_gens[nscat]->getPosition() ;
         winv = 1.0/target_gens[nscat]->getWeight(pos) ;
         return pos ;
       } else {
         // scattering vertex
         return generateScatteringVertexPosition(nscat,nv,winv) ;
       }
     }

◆ generateScatteringVertexPosition()

virtual JPosition3D JMARKOV::JMarkovEnsembleIntegrator::generateScatteringVertexPosition	(	int	nscat,
		int	nv,
		double &	winv
	)

protectedpure virtual

Implemented in JMARKOV::JMarkovEnsembleIntegrator3D, and JMARKOV::JMarkovEnsembleIntegrator1D.

◆ getIndex()

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

inlineprotected

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

Definition at line 251 of file JMarkovIntegrator.hh.

                                       {
       // indices are distributed as follows
       // [0]: nscat=1, nv=1
       // [1]: nscat=2, nv=1
       // [2]: nscat=2, nv=2
       // [3]: nscat=3, nv=1
       // [4]: nscat=3, nv=2
       // etc.
       return (nscat*(nscat-1))/2 + nv - 1 ;
     }

◆ initgenerators()

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

inlineprotected

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.

                                      {
       // target vertex
       init_target_generator(nscat) ;
  
       // scattering vertices
       init_scattering_vertex_generators(nscat) ;
     }

◆ init_target_generator()

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

inlineprotected

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

Definition at line 272 of file JMarkovIntegrator.hh.

                                             {
       if( (int)target_gens.size() <= nscat ) target_gens.resize( nscat+1 ) ;
  
       if( trg->getRadius()>0 ) {
         // for a finite target, we fill a "sphere generator" from the ensemble
         char hname[200] ;
         sprintf( hname, "htarget_nscat%i", nscat ) ;
         TH2D* ht = new TH2D(hname,";cos(#theta);#phi",100,-1,1,100,-M_PI,M_PI) ;
         for( vector<JPhotonPath>::iterator it=ensemble.begin() ; it!=ensemble.end() ; ++it ) {
           if( it->n-2 != nscat ) continue ;
           JVersor3D dir( it->back() - trg->getPosition() ) ;
           ht->Fill( dir.getDZ(), dir.getPhi() ) ;
         }
         target_gens[nscat] = new JSphereGenerator( trg->getPosition(), trg->getRadius(), ht ) ;
         delete ht ;
       } else {
         // for an infinitesimal target, we create a "generator" that only generates the single point
         target_gens[nscat] = new JSphereGenerator( trg->getPosition() ) ;
       }
     }

◆ init_scattering_vertex_generators()

virtual void JMARKOV::JMarkovEnsembleIntegrator::init_scattering_vertex_generators ( int nscat )

protectedpure virtual

Implemented in JMARKOV::JMarkovEnsembleIntegrator3D, and JMARKOV::JMarkovEnsembleIntegrator1D.

◆ free_target_gens()

void JMARKOV::JMarkovEnsembleIntegrator::free_target_gens ( )

inlineprotected

Definition at line 295 of file JMarkovIntegrator.hh.

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

◆ 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.

                                                                                                                                                {
     vector<double> contributions(N,-1) ;
  
     for( int i=0 ; i<N ; ++i ) {
       double winv ;
       JPhotonPath p = generatePath(nscat,winv) ;
  
       double rho = getPhotonPathProbabilityDensity(p,src,sm,trg,lambda_abs) ;
       contributions[i] = rho * winv ; 
     }
     return contributions ;
   }

◆ 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.

                                                                                                                                                      {
     vector<double> contributions(N,-1) ;
  
     for( int i=0 ; i<N ; ++i ) {
       const double r = 10 ;
       double winv ;
       JPhotonPath p = generatePath(nscat,winv) ;
       double rho = 1.0/(4.0/3.0*M_PI*r*r*r) ;
       if( p[1].getLength()>r ) rho = 0 ;
       contributions[i] = rho * winv ; 
     }
     return contributions ;
   }

◆ 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.

                                                                                                                                                                   {
     vector<JPhotonPath> paths ;
  
     for( int i=0 ; i<N ; ++i ) {
       double winv ;
       JPhotonPath p = generatePath(nscat,winv) ;
       paths.push_back(p) ;
     }
     return paths ;
   }

◆ 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.

                                                                 {
       // this default implementation assumes that the vertex positions
       // are completely uncorrelated
       JPhotonPath p(nscat) ;
       double _winv = 1 ;
       for( int nv=0 ; nv<nscat+2 ; ++nv ) {
         double part_winv ;
         p[nv] = generatePosition( nscat, nv, part_winv ) ;
         _winv *= part_winv ;
       }
       winv = _winv ;
       return p ;
     }

Member Data Documentation

◆ target_gens

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

protected

Definition at line 301 of file JMarkovIntegrator.hh.

◆ ensemble

vector<JPhotonPath> JMARKOV::JMarkovEnsembleIntegrator::ensemble

protected

Definition at line 302 of file JMarkovIntegrator.hh.

◆ nbinsx

int JMARKOV::JMarkovEnsembleIntegrator::nbinsx

protected

Definition at line 303 of file JMarkovIntegrator.hh.

◆ nbinsy

int JMARKOV::JMarkovEnsembleIntegrator::nbinsy

protected

Definition at line 304 of file JMarkovIntegrator.hh.

◆ nbinsz

int JMARKOV::JMarkovEnsembleIntegrator::nbinsz

protected

Definition at line 305 of file JMarkovIntegrator.hh.

◆ trg

JTargetModel* JMARKOV::JMarkovEnsembleIntegrator::trg

protected

Definition at line 306 of file JMarkovIntegrator.hh.

◆ source_position

JPosition3D JMARKOV::JMarkovEnsembleIntegrator::source_position

protected

Definition at line 307 of file JMarkovIntegrator.hh.

The documentation for this class was generated from the following file:

JMarkovIntegrator.hh

Public Member Functions

Protected Member Functions

Protected Attributes

Detailed Description

Constructor & Destructor Documentation

◆ JMarkovEnsembleIntegrator()

◆ ~JMarkovEnsembleIntegrator()

Member Function Documentation

◆ writeHistograms()

◆ generatePosition()

◆ generateScatteringVertexPosition()

◆ getIndex()

◆ initgenerators()

◆ init_target_generator()

◆ init_scattering_vertex_generators()

◆ free_target_gens()

◆ integrate()

◆ dummy_integrate()

◆ get_diagnostic_ensemble()

◆ generatePath()

Member Data Documentation

◆ target_gens

◆ ensemble

◆ nbinsx

◆ nbinsy

◆ nbinsz

◆ trg

◆ source_position