G4ITModelProcessor Class Reference

#include <G4ITModelProcessor.hh>


Public Member Functions

 G4ITModelProcessor ()
virtual ~G4ITModelProcessor ()
void SetModelHandler (G4ITModelHandler *)
void Initialize ()
void CleanProcessor ()
void InitializeStepper (const G4double &currentGlobalTime, const G4double &userMinTime)
void CalculateTimeStep (const G4Track *, const G4double)
void DoCalculateStep ()
void FindReaction (std::map< G4Track *, G4TrackVectorHandle > *, const double currentStepTime, const double previousStepTime, const bool reachedUserStepTimeLimit)
const std::vector< std::vector<
G4VITModel * > > * 
GetCurrentModel ()
std::vector< G4ITReactionChange * > * GetReactionInfo ()
const G4TrackGetTrack () const

Protected Member Functions

void SetTrack (const G4Track *)
 G4ITModelProcessor (const G4ITModelProcessor &other)
G4ITModelProcessoroperator= (const G4ITModelProcessor &other)

Protected Attributes

G4bool fInitialized
G4ITModelHandlerfpModelHandler
const G4TrackfpTrack
G4double fUserMinTimeStep
std::vector< std::vector<
G4VITModel * > > 
fCurrentModel
G4VITModelfpModel
G4ITModelManagerfpModelManager
G4ITType fCurrentType1
G4ITType fCurrentType2
std::vector< G4ITReactionChange * > fReactionInfo

Static Protected Attributes

static std::map< const G4Track *,
G4bool
fHasReacted


Detailed Description

The G4ITModelProcessor will call the two processes defined in G4VITModel. This processes act at the beginning and end of each step. The first one, the TimeStepper will calculate a time step to propagate all the track and eventually it can return some tracks that can likely react at the end of the step. The second one, the ReactionProcess will make the tracks reacting.

Definition at line 62 of file G4ITModelProcessor.hh.


Constructor & Destructor Documentation

G4ITModelProcessor::G4ITModelProcessor (  ) 

Default constructor

Definition at line 42 of file G4ITModelProcessor.cc.

References fCurrentModel, fInitialized, fpModel, fpModelHandler, fpModelManager, fpTrack, fUserMinTimeStep, and G4ITType::size().

00043 {
00044     //ctor
00045     fpTrack = 0;
00046     fpModelHandler = 0;
00047     fpModel = 0;
00048     fInitialized = false;
00049     fpModelManager = 0;
00050     fCurrentModel.assign(G4ITType::size(), std::vector<G4VITModel*>());
00051 
00052     for(int i = 0 ; i < (int) G4ITType::size() ; i++)
00053     {
00054         fCurrentModel[i].assign(G4ITType::size(),0);
00055     }
00056     fUserMinTimeStep = -1.;
00057 }

G4ITModelProcessor::~G4ITModelProcessor (  )  [virtual]

Default destructor

Definition at line 59 of file G4ITModelProcessor.cc.

References fCurrentModel, and fReactionInfo.

00060 {
00061     //dtor
00062 //    if(fpModelHandler) delete fpModelHandler; deleted by G4ITStepManager
00063     fCurrentModel.clear();
00064     fReactionInfo.clear();
00065 }

G4ITModelProcessor::G4ITModelProcessor ( const G4ITModelProcessor other  )  [protected]

Copy constructor

Parameters:
other Object to copy from

Definition at line 68 of file G4ITModelProcessor.cc.

References fInitialized, fpModel, fpModelHandler, fpModelManager, fpTrack, and fUserMinTimeStep.

00069 {
00070     //copy ctorr
00071     fpTrack = 0;
00072     fpModelHandler = 0;
00073     fpModel = 0;
00074     fInitialized = false;
00075     fpModelManager = 0;
00076     fUserMinTimeStep = -1.;
00077 }


Member Function Documentation

void G4ITModelProcessor::CalculateTimeStep ( const G4Track ,
const   G4double 
)

Definition at line 156 of file G4ITModelProcessor.cc.

References CleanProcessor(), DoCalculateStep(), FatalErrorInArgument, fUserMinTimeStep, G4Exception(), and SetTrack().

00157 {
00158     // G4cout  << "G4ITModelProcessor::CalculateStep" << G4endl;
00159     CleanProcessor();
00160     if(track == 0)
00161     {
00162         G4ExceptionDescription exceptionDescription ;
00163         exceptionDescription << "No track was passed to the method (track == 0).";
00164         G4Exception("G4ITModelProcessor::CalculateStep","ITModelProcessor004",
00165                     FatalErrorInArgument,exceptionDescription);
00166     }
00167     SetTrack(track);
00168     fUserMinTimeStep = userMinTimeStep ;
00169 
00170     DoCalculateStep();
00171 }

void G4ITModelProcessor::CleanProcessor (  )  [inline]

Restaure original state of the modelProcessor. This method should be call only by the ITStepManager

Definition at line 178 of file G4ITModelProcessor.hh.

References fpTrack.

Referenced by CalculateTimeStep().

00179 {
00180     fpTrack = 0;
00181 }

void G4ITModelProcessor::DoCalculateStep (  ) 

Definition at line 174 of file G4ITModelProcessor.cc.

References fCurrentModel, fpModel, fpTrack, fUserMinTimeStep, GetIT(), and G4IT::GetITType().

Referenced by CalculateTimeStep().

00175 {
00176     if(fpModel) // ie only one model has been declared and will be used
00177     {
00178         fpModel -> GetTimeStepper()->CalculateStep(*fpTrack, fUserMinTimeStep);
00179     }
00180     else // ie many models have been declared and will be used
00181     {
00182         std::vector<G4VITModel*>& model = fCurrentModel[GetIT(fpTrack)->GetITType()];
00183 
00184         for(int i =0 ; i < (int) model.size() ; i++)
00185         {
00186             if(model[i] == 0) continue;
00187             model[i]->GetTimeStepper()->CalculateStep(*fpTrack, fUserMinTimeStep);
00188         }
00189     }
00190 }

void G4ITModelProcessor::FindReaction ( std::map< G4Track *, G4TrackVectorHandle > *  ,
const double  currentStepTime,
const double  previousStepTime,
const bool  reachedUserStepTimeLimit 
)

Definition at line 193 of file G4ITModelProcessor.cc.

References FatalErrorInArgument, fCurrentModel, fHasReacted, fpModelHandler, fReactionInfo, fStopAndKill, G4Exception(), G4ITModelHandler::GetAllModelManager(), GetIT(), G4Track::GetTrackStatus(), G4VITReactionProcess::MakeReaction(), and G4VITReactionProcess::ResetChanges().

00197 {
00198     // DEBUG
00199     //    G4cout << "G4ITReactionManager::FindReaction" << G4endl;
00200     if(tracks == 0)       return ;
00201 
00202     if(fpModelHandler->GetAllModelManager()->empty()) return ;
00203 
00204     std::map<G4Track*, G4TrackVectorHandle>::iterator tracks_i = tracks->begin();;
00205 
00206     for(tracks_i = tracks->begin() ; tracks_i != tracks-> end() ; tracks_i ++)
00207     {
00209         G4Track* trackA = tracks_i->first;
00210 
00211         if(trackA == 0)         continue;
00212 
00213         std::map<const G4Track*, G4bool>::iterator it_hasReacted = fHasReacted.find(trackA);
00214         if(it_hasReacted != fHasReacted.end()) continue;
00215         if(trackA->GetTrackStatus() == fStopAndKill) continue;
00216 
00217         G4IT* ITA = GetIT(trackA);
00218         G4ITType ITypeA = ITA -> GetITType();
00219 
00220         const std::vector<G4VITModel*> model = fCurrentModel[ITypeA];
00221 
00222         G4TrackVectorHandle& trackB_vector = tracks_i->second ;
00223         std::vector<G4Track*>::iterator trackB_i = trackB_vector->begin();
00224 
00225         G4Track* trackB = 0 ;
00226         G4ITType ITypeB(-1);
00227         G4VITReactionProcess* process = 0;
00228         G4ITReactionChange* changes = 0;
00229 
00230         for(; trackB_i != trackB_vector->end() ; trackB_i++)
00231         {
00232             trackB = *trackB_i;
00233 
00234             if(trackB == 0)         continue;
00235             it_hasReacted = fHasReacted.find(trackB);
00236             if(it_hasReacted != fHasReacted.end()) continue;
00237             if(trackB->GetTrackStatus() == fStopAndKill) continue;
00238 
00239             // DEBUG
00240             //             G4cout << "Couple : " << trackA->GetParticleDefinition->GetParticleName() << " ("
00241             //                        << trackA->GetTrackID() << ")   "
00242             //                        << trackB->GetParticleDefinition->GetParticleName() << " ("
00243             //                        << trackB->GetTrackID() << ")"
00244             //                        << G4endl;
00245 
00246             if(trackB == trackA)
00247             {
00248                 G4ExceptionDescription exceptionDescription ;
00249                 exceptionDescription << "The IT reaction process sent back a reaction between trackA and trackB. ";
00250                 exceptionDescription << "The problem is trackA == trackB";
00251                 G4Exception("G4ITModelProcessor::FindReaction","ITModelProcessor005",
00252                             FatalErrorInArgument,exceptionDescription);
00253             }
00254 
00255             G4IT* ITB = GetIT(trackB);
00256             G4ITType ITypeBtmp = ITB -> GetITType();
00257 
00258             if(ITypeB != ITypeBtmp)
00259             {
00260                 ITypeB = ITypeBtmp ;
00261 
00262                 if(model[ITypeB])
00263                     process = model[ITypeB]->GetReactionProcess();
00264             }
00265 
00266             if(process && process -> TestReactibility(*trackA, *trackB,
00267                                                       currentStepTime, previousStepTime,
00268                                                       reachedUserStepTimeLimit))
00269             {
00270                 changes = process->MakeReaction(*trackA, *trackB);
00271             }
00272 
00273             if(changes)
00274             {
00275                 fHasReacted[trackA] = true;
00276                 fHasReacted[trackB] = true;
00277                 changes -> GetTrackA();
00278                 changes -> GetTrackB();
00279 
00280                 fReactionInfo.push_back(changes);
00281 
00282                 process->ResetChanges();
00283                 changes = 0;
00284 
00285                 break;
00286             }
00287         }
00288     }
00289 
00290     fHasReacted.clear();
00291 }

const std::vector< std::vector< G4VITModel * > > * G4ITModelProcessor::GetCurrentModel (  )  [inline]

Definition at line 161 of file G4ITModelProcessor.hh.

References fCurrentModel.

00162 {
00163     return &fCurrentModel ;
00164 }

std::vector<G4ITReactionChange*>* G4ITModelProcessor::GetReactionInfo (  )  [inline]

Definition at line 104 of file G4ITModelProcessor.hh.

References fReactionInfo.

00105     {
00106         return &fReactionInfo;
00107     }

const G4Track* G4ITModelProcessor::GetTrack (  )  const [inline]

Definition at line 109 of file G4ITModelProcessor.hh.

References fpTrack.

00110     {
00111         return fpTrack;
00112     }

void G4ITModelProcessor::Initialize (  ) 

Definition at line 87 of file G4ITModelProcessor.cc.

References fInitialized, fpModelHandler, and G4ITModelHandler::Initialize().

00088 {
00089     fpModelHandler->Initialize();
00090     fInitialized = true;
00091 }

void G4ITModelProcessor::InitializeStepper ( const G4double currentGlobalTime,
const G4double userMinTime 
)

Definition at line 94 of file G4ITModelProcessor.cc.

References FatalErrorInArgument, fCurrentModel, fpModel, fpModelHandler, fpModelManager, G4Exception(), G4ITModelHandler::GetAllModelManager(), G4VITModel::GetTimeStepper(), and G4VITTimeStepper::SetTimes().

00096 {
00097     // G4cout << "G4ITModelProcessor::InitializeStepper" << G4endl;
00098     if(fpModelHandler==0)
00099     {
00100         G4ExceptionDescription exceptionDescription ;
00101         exceptionDescription << "No G4ITModelHandler was passed to the modelProcessor.";
00102         G4Exception("G4ITModelProcessor::InitializeStepper","ITModelProcessor002",
00103                     FatalErrorInArgument,exceptionDescription);
00104     }
00105     const std::vector<std::vector<G4ITModelManager*> >* modelManager = fpModelHandler
00106             ->GetAllModelManager();
00107 
00108     if(modelManager==0)
00109     {
00110         G4ExceptionDescription exceptionDescription ;
00111         exceptionDescription << "No G4ITModelManager was register to G4ITModelHandler.";
00112         G4Exception("G4ITModelProcessor::InitializeStepper","ITModelProcessor003",
00113                     FatalErrorInArgument,exceptionDescription);
00114     }
00115 
00116     int nbModels1 = modelManager->size() ;
00117 
00118     G4VITTimeStepper::SetTimes(currentGlobalTime, userMinTime) ;
00119 
00120     // TODO !!!
00121     //    if( nbModels1 != 1 || (nbModels1 == 1 && !fpModelManager) )
00122     {
00123         int nbModels2 = -1;
00124         G4VITModel* model = 0;
00125         G4ITModelManager* modman = 0;
00126 
00127         for(int i = 0 ; i < nbModels1 ; i++)
00128         {
00129             nbModels2 = (*modelManager)[i].size();
00130 
00131             for(int j = 0 ; j <= i ; j++)
00132             {
00133                 modman = (*modelManager)[i][j];
00134 
00135                 if(modman == 0) continue ;
00136 
00137                 model       =  modman -> GetModel(currentGlobalTime);
00138                 G4VITTimeStepper* stepper   = model->GetTimeStepper() ;
00139 
00140 //                stepper -> PrepareForAllProcessors() ;
00141                 stepper -> Prepare() ;
00142                 fCurrentModel[i][j] = model;
00143             }
00144         }
00145 
00146         if(nbModels1 == 1 && nbModels2 ==1)
00147         {
00148             fpModelManager = modman;
00149             fpModel = model;
00150         }
00151         else fpModel = 0;
00152     }
00153 }

G4ITModelProcessor & G4ITModelProcessor::operator= ( const G4ITModelProcessor other  )  [protected]

Assignment operator

Parameters:
other Object to assign from
Returns:
A reference to this

Definition at line 80 of file G4ITModelProcessor.cc.

00081 {
00082     if (this == &rhs) return *this; // handle self assignment
00083     //assignment operator
00084     return *this;
00085 }

void G4ITModelProcessor::SetModelHandler ( G4ITModelHandler  )  [inline]

Definition at line 166 of file G4ITModelProcessor.hh.

References FatalErrorInArgument, fInitialized, fpModelHandler, and G4Exception().

00167 {
00168     if(fInitialized == 1)
00169     {
00170         G4ExceptionDescription exceptionDescription ;
00171         exceptionDescription << "You are trying to set a new model while the model processor has alreaday be initialized";
00172         G4Exception("G4ITModelProcessor::SetModelHandler","ITModelProcessor001",
00173                     FatalErrorInArgument,exceptionDescription);
00174     }
00175     fpModelHandler = modelHandler;
00176 }

void G4ITModelProcessor::SetTrack ( const G4Track  )  [inline, protected]

Definition at line 156 of file G4ITModelProcessor.hh.

References fpTrack.

Referenced by CalculateTimeStep().

00157 {
00158     fpTrack = track;
00159 }


Field Documentation

std::vector<std::vector<G4VITModel*> > G4ITModelProcessor::fCurrentModel [protected]

Definition at line 136 of file G4ITModelProcessor.hh.

Referenced by DoCalculateStep(), FindReaction(), G4ITModelProcessor(), GetCurrentModel(), InitializeStepper(), and ~G4ITModelProcessor().

G4ITType G4ITModelProcessor::fCurrentType1 [protected]

Definition at line 144 of file G4ITModelProcessor.hh.

G4ITType G4ITModelProcessor::fCurrentType2 [protected]

Definition at line 145 of file G4ITModelProcessor.hh.

std::map< const G4Track *, G4bool > G4ITModelProcessor::fHasReacted [static, protected]

Definition at line 149 of file G4ITModelProcessor.hh.

Referenced by FindReaction().

G4bool G4ITModelProcessor::fInitialized [protected]

Definition at line 128 of file G4ITModelProcessor.hh.

Referenced by G4ITModelProcessor(), Initialize(), and SetModelHandler().

G4VITModel* G4ITModelProcessor::fpModel [protected]

Definition at line 140 of file G4ITModelProcessor.hh.

Referenced by DoCalculateStep(), G4ITModelProcessor(), and InitializeStepper().

G4ITModelHandler* G4ITModelProcessor::fpModelHandler [protected]

Definition at line 129 of file G4ITModelProcessor.hh.

Referenced by FindReaction(), G4ITModelProcessor(), Initialize(), InitializeStepper(), and SetModelHandler().

G4ITModelManager* G4ITModelProcessor::fpModelManager [protected]

Definition at line 141 of file G4ITModelProcessor.hh.

Referenced by G4ITModelProcessor(), and InitializeStepper().

const G4Track* G4ITModelProcessor::fpTrack [protected]

Definition at line 131 of file G4ITModelProcessor.hh.

Referenced by CleanProcessor(), DoCalculateStep(), G4ITModelProcessor(), GetTrack(), and SetTrack().

std::vector<G4ITReactionChange*> G4ITModelProcessor::fReactionInfo [protected]

Definition at line 148 of file G4ITModelProcessor.hh.

Referenced by FindReaction(), GetReactionInfo(), and ~G4ITModelProcessor().

G4double G4ITModelProcessor::fUserMinTimeStep [protected]

Definition at line 132 of file G4ITModelProcessor.hh.

Referenced by CalculateTimeStep(), DoCalculateStep(), and G4ITModelProcessor().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:19 2013 for Geant4 by  doxygen 1.4.7