445 lines
14 KiB
C++
445 lines
14 KiB
C++
// Written by Peter Stromberger
|
|
// based on work of Mirco Deckenhoff
|
|
// modified by Bastian Rössler
|
|
|
|
#include "Analysis.hh"
|
|
#include "Parameters.hh"
|
|
#include <ctime>
|
|
#include <fstream>
|
|
|
|
|
|
Analysis* Analysis::singleton = 0;
|
|
|
|
Analysis::Analysis()
|
|
{
|
|
// Set file name for simulation output
|
|
G4int threadNumber;
|
|
std::ifstream inFile;
|
|
inFile.open("threadNumber.txt");
|
|
inFile >> threadNumber;
|
|
inFile.close();
|
|
|
|
sprintf(fileName,"outFile_%i.root", threadNumber);
|
|
|
|
detectedPhotons = new TTree("DetectedPhotons","Photons detected at the detector strip.");
|
|
|
|
detectedPhotons->Branch("runID", &runIDBuffer, "runID/F"); //modifed by me
|
|
detectedPhotons->Branch("eventID", &eventIDBuffer, "eventID/F"); //modifed by me
|
|
detectedPhotons->Branch("detNumb", &detNumbBuffer, "detNumb/F"); //modifed by me
|
|
|
|
detectedPhotons->Branch("xPixel", &xPixelBuffer, "xPixel/F"); //modifed by me
|
|
detectedPhotons->Branch("yPixel", &yPixelBuffer, "yPixel/F"); //modifed by me
|
|
detectedPhotons->Branch("energy",&energyBuffer,"energy/F");
|
|
detectedPhotons->Branch("wavelength",&wavelengthBuffer,"wavelength/F");
|
|
detectedPhotons->Branch("localtime",&timeBuffer,"localtime/F");
|
|
detectedPhotons->Branch("abstime",&absTimeBuffer,"abstime/F");
|
|
detectedPhotons->Branch("length",&lengthBuffer,"length/F");
|
|
detectedPhotons->Branch("x",&xBuffer,"x/F");
|
|
detectedPhotons->Branch("y",&yBuffer,"y/F");
|
|
detectedPhotons->Branch("z",&zBuffer,"z/F");
|
|
detectedPhotons->Branch("px",&pxBuffer,"px/F");
|
|
detectedPhotons->Branch("py",&pyBuffer,"py/F");
|
|
detectedPhotons->Branch("pz",&pzBuffer,"pz/F");
|
|
detectedPhotons->Branch("vertexX",&vertexXBuffer,"vertexX/F");
|
|
detectedPhotons->Branch("vertexY",&vertexYBuffer,"vertexY/F");
|
|
detectedPhotons->Branch("vertexZ",&vertexZBuffer,"vertexZ/F");
|
|
detectedPhotons->Branch("vertexPx",&vertexPxBuffer,"vertexPx/F");
|
|
detectedPhotons->Branch("vertexPy",&vertexPyBuffer,"vertexPy/F");
|
|
detectedPhotons->Branch("vertexPz",&vertexPzBuffer,"vertexPz/F");
|
|
detectedPhotons->Branch("gpsPosX",&gpsPositionX,"gpsPosX/F");
|
|
detectedPhotons->Branch("gpsPosY",&gpsPositionY,"gpsPosY/F");
|
|
detectedPhotons->Branch("gpsPosZ",&gpsPositionZ,"gpsPosZ/F");
|
|
detectedPhotons->Branch("gpsDirX",&gpsDirectionX,"gpsPosX/F");
|
|
detectedPhotons->Branch("gpsDirY",&gpsDirectionY,"gpsPosY/F");
|
|
detectedPhotons->Branch("gpsDirZ",&gpsDirectionZ,"gpsPosZ/F");
|
|
detectedPhotons->Branch("runId",&runIdBuffer,"runId/I");
|
|
detectedPhotons->Branch("eventId",&eventIdBuffer,"eventId/I");
|
|
detectedPhotons->Branch("trackId",&trackIdBuffer,"trackId/I");
|
|
detectedPhotons->Branch("creatorProcess",&creatorProcessBuffer,"creatorProcess/I");
|
|
detectedPhotons->Branch("parentId",&parentIdBuffer,"parentId/I");
|
|
detectedPhotons->Branch("reflMirr",&reflMirrBuffer,"reflMirr/I");
|
|
detectedPhotons->Branch("reflSurf",&reflSurfBuffer,"reflSurf/I");
|
|
detectedPhotons->Branch("reflTotalCladClad",&reflTotalCladCladBuffer,"reflTotalCladClad/I");
|
|
detectedPhotons->Branch("reflTotalCoreClad",&reflTotalCoreCladBuffer,"reflTotalCoreClad/I");
|
|
detectedPhotons->Branch("reflFresnelCladClad",&reflFresnelCladCladBuffer,"reflFresnelCladClad/I");
|
|
detectedPhotons->Branch("reflFresnelCoreClad",&reflFresnelCoreCladBuffer,"reflFresnelCoreClad/I");
|
|
detectedPhotons->Branch("refracCladClad",&refracCladCladBuffer,"refracCladClad/I");
|
|
detectedPhotons->Branch("refracCoreClad",&refracCoreCladBuffer,"refracCoreClad/I");
|
|
detectedPhotons->Branch("rayleighScatterings",&rayleighScatteringsBuffer,"rayleighScatterings/I");
|
|
detectedPhotons->Branch("lengthInCore",&lengthInCoreBuffer,"lengthInCore/F");
|
|
detectedPhotons->Branch("lengthInInnerCladding",&lengthInInnerCladdingBuffer,"lengthInInnerCladding/F");
|
|
detectedPhotons->Branch("lengthInOuterCladding",&lengthInOuterCladdingBuffer,"lengthInOuterCladding/F");
|
|
|
|
trigger = new TTree("Trigger","Hits in the trigger.");
|
|
|
|
trigger->Branch("runID", &runIDBuffer, "runID/F");
|
|
trigger->Branch("eventID", &eventIDBuffer, "eventID/F");
|
|
trigger->Branch("edep", &edepBuffer, "edep/F");
|
|
trigger->Branch("xPos", &xPosBuffer, "xPos/F");
|
|
trigger->Branch("yPos", &yPosBuffer, "yPos/F");
|
|
trigger->Branch("zPos", &zPosBuffer, "zPos/F");
|
|
|
|
energyTrack = new TTree("EnergyTrack", "Deposited energy and tracklength in core of fibre.");
|
|
|
|
energyTrack->Branch("runID", &runIDBuffer, "runID/F");
|
|
energyTrack->Branch("eventID", &eventIDBuffer, "eventID/F");
|
|
energyTrack->Branch("edep", &energyBuffer, "edep/F");
|
|
energyTrack->Branch("trackL", &lengthBuffer, "trackL/F");
|
|
|
|
h_energy = new TH1F("h_energy","h_energy",200,0.0,1.0);
|
|
h_trackL = new TH1F("h_trackL","h_trackL",400,0.0,2.0);
|
|
|
|
primaryParticleTrack = new TTree("PrimaryParticleTrack", "Track of the primary particles");
|
|
|
|
primaryParticleTrack->Branch("runID", &runIDBuffer, "runID/F");
|
|
primaryParticleTrack->Branch("eventID", &eventIDBuffer, "eventID/F");
|
|
primaryParticleTrack->Branch("xPos", &xPosBuffer, "xPos/F");
|
|
primaryParticleTrack->Branch("yPos", &yPosBuffer, "yPos/F");
|
|
primaryParticleTrack->Branch("zPos", &zPosBuffer, "zPos/F");
|
|
|
|
initialParticle = new TTree("InitialParticle", "Initial Particle");
|
|
|
|
initialParticle->Branch("runID", &runIDBuffer, "runID/F");
|
|
initialParticle->Branch("eventID", &eventIDBuffer, "eventID/F");
|
|
initialParticle->Branch("energy", &energyBuffer, "energy/F");
|
|
initialParticle->Branch("xMom", &xMomBuffer, "xMom/F");
|
|
initialParticle->Branch("yMom", &yMomBuffer, "yMom/F");
|
|
initialParticle->Branch("zMom", &zMomBuffer, "zMom/F");
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
void Analysis::PrepareNewRun(const G4Run* aRun)
|
|
{
|
|
runIdBuffer = aRun->GetRunID();
|
|
}
|
|
|
|
|
|
void Analysis::EndOfRun()
|
|
{
|
|
dataFile = new TFile(fileName,"update");
|
|
|
|
//detectedPhotons->Write("",TObject::kOverwrite);
|
|
//trigger->Write("", TObject::kOverwrite);
|
|
// primaryParticleTrack->Write("", TObject::kOverwrite);
|
|
//initialParticle->Write("", TObject::kOverwrite);
|
|
// energyTrack->Write("", TObject::kOverwrite);
|
|
h_energy->Write("", TObject::kOverwrite);
|
|
h_trackL->Write("", TObject::kOverwrite);
|
|
G4cout << "Data written to file " << fileName << G4endl;
|
|
|
|
dataFile->Close();
|
|
|
|
G4cout << "Data file closed: " << fileName << G4endl;
|
|
|
|
delete dataFile;
|
|
}
|
|
|
|
|
|
void Analysis::Close()
|
|
{
|
|
delete detectedPhotons;
|
|
delete trigger;
|
|
delete initialParticle;
|
|
delete energyTrack;
|
|
}
|
|
|
|
|
|
void Analysis::PrepareNewEvent(const G4Event* anEvent)
|
|
{
|
|
eventIdBuffer = anEvent->GetEventID();
|
|
|
|
reflectionsAtMirror.resize(0);
|
|
reflectionsAtFibreSurface.resize(0);
|
|
reflectionsTotalAtCladCladInterface.resize(0);
|
|
reflectionsTotalAtCoreCladInterface.resize(0);
|
|
reflectionsFresnelAtCladCladInterface.resize(0);
|
|
reflectionsFresnelAtCoreCladInterface.resize(0);
|
|
|
|
refractionsAtCladCladInterface.resize(0);
|
|
refractionsAtCoreCladInterface.resize(0);
|
|
|
|
rayleighScatterings.resize(0);
|
|
|
|
lengthInCore.resize(0);
|
|
lengthInInnerCladding.resize(0);
|
|
lengthInOuterCladding.resize(0);
|
|
}
|
|
|
|
void Analysis::FillInitialParticle(Double_t runID, Double_t eventID, Double_t energy, Double_t xMom, Double_t yMom, Double_t zMom)
|
|
{
|
|
if(initialParticle != NULL)
|
|
{
|
|
runIDBuffer = runID;
|
|
eventIDBuffer = eventID;
|
|
energyBuffer = energy;
|
|
xMomBuffer = xMom;
|
|
yMomBuffer = yMom;
|
|
zMomBuffer = zMom;
|
|
// initialParticle->Fill();
|
|
}
|
|
else
|
|
G4cout << "Nullpointer to initialParticle tree." << G4endl;
|
|
}
|
|
|
|
void Analysis::FillTrigger(Double_t runID, Double_t eventID, Double_t edep, Double_t xPos, Double_t yPos, Double_t zPos)
|
|
{
|
|
if(trigger != NULL)
|
|
{
|
|
runIDBuffer = runID;
|
|
eventIDBuffer = eventID;
|
|
edepBuffer = edep;
|
|
xPosBuffer = xPos;
|
|
yPosBuffer = yPos;
|
|
zPosBuffer = zPos;
|
|
|
|
// trigger->Fill();
|
|
}
|
|
else
|
|
{
|
|
G4cout << "Nullpointer to trigger tree." << G4endl;
|
|
}
|
|
}
|
|
|
|
void Analysis::FillEnergyTrack(Double_t runID, Double_t eventID, Double_t energy, Double_t trackL)
|
|
{
|
|
if(energyTrack != NULL)
|
|
{
|
|
runIDBuffer = runID;
|
|
eventIDBuffer = eventID;
|
|
energyBuffer = energy;
|
|
lengthBuffer = trackL;
|
|
// energyTrack->Fill();
|
|
h_energy->Fill(energy);
|
|
h_trackL->Fill(trackL);
|
|
}
|
|
else
|
|
{
|
|
G4cout << "Nullpointer to energyTrack tree" << G4endl;
|
|
}
|
|
}
|
|
|
|
void Analysis::FillPrimaryParticleTrack(Double_t runID, Double_t eventID,
|
|
Double_t xPos, Double_t yPos, Double_t zPos)
|
|
{
|
|
if(primaryParticleTrack != NULL)
|
|
{
|
|
runIDBuffer = runID;
|
|
eventIDBuffer = eventID;
|
|
xPosBuffer = xPos;
|
|
yPosBuffer = yPos;
|
|
zPosBuffer = zPos;
|
|
// primaryParticleTrack->Fill();
|
|
}
|
|
else
|
|
{
|
|
G4cout << "Nullpointer to primaryParticleTrack tree" << G4endl;
|
|
}
|
|
}
|
|
|
|
void Analysis::FillDetectedPhotons(Double_t runID, Double_t eventID, Double_t detNumb, Double_t xPixel,
|
|
Double_t yPixel, Double_t energy,
|
|
Double_t time, Double_t length, Double_t absTime,
|
|
Double_t x, Double_t y, Double_t z,
|
|
Double_t px, Double_t py, Double_t pz,
|
|
Double_t vertexX, Double_t vertexY, Double_t vertexZ,
|
|
Double_t vertexPx, Double_t vertexPy, Double_t vertexPz,
|
|
Int_t trackId,
|
|
Int_t creatorProcess, Int_t parentId)
|
|
{
|
|
if(detectedPhotons!=NULL)
|
|
{
|
|
runIDBuffer = runID;
|
|
eventIDBuffer = eventID;
|
|
detNumbBuffer = detNumb;
|
|
xPixelBuffer = xPixel;
|
|
yPixelBuffer = yPixel;
|
|
|
|
energyBuffer = energy;
|
|
timeBuffer = time;
|
|
lengthBuffer = length;
|
|
|
|
wavelengthBuffer = 1239.842/energy;
|
|
|
|
absTimeBuffer = absTime;
|
|
xBuffer = x;
|
|
yBuffer = y;
|
|
zBuffer = z;
|
|
|
|
pxBuffer = px;
|
|
pyBuffer = py;
|
|
pzBuffer = pz;
|
|
|
|
vertexXBuffer = vertexX;
|
|
vertexYBuffer = vertexY;
|
|
vertexZBuffer = vertexZ;
|
|
|
|
vertexPxBuffer = vertexPx;
|
|
vertexPyBuffer = vertexPy;
|
|
vertexPzBuffer = vertexPz;
|
|
|
|
trackIdBuffer = trackId;
|
|
|
|
creatorProcessBuffer = creatorProcess;
|
|
parentIdBuffer = parentId;
|
|
|
|
reflMirrBuffer = reflectionsAtMirror[trackId];
|
|
reflSurfBuffer = reflectionsAtFibreSurface[trackId];
|
|
reflTotalCladCladBuffer = reflectionsTotalAtCladCladInterface[trackId];
|
|
reflTotalCoreCladBuffer = reflectionsTotalAtCoreCladInterface[trackId];
|
|
reflFresnelCladCladBuffer = reflectionsFresnelAtCladCladInterface[trackId];
|
|
reflFresnelCoreCladBuffer = reflectionsFresnelAtCoreCladInterface[trackId];
|
|
|
|
refracCladCladBuffer = refractionsAtCladCladInterface[trackId];
|
|
refracCoreCladBuffer = refractionsAtCoreCladInterface[trackId];
|
|
|
|
rayleighScatteringsBuffer = rayleighScatterings[trackId];
|
|
|
|
lengthInCoreBuffer = lengthInCore[trackId];
|
|
lengthInInnerCladdingBuffer = lengthInInnerCladding[trackId];
|
|
lengthInOuterCladdingBuffer = lengthInOuterCladding[trackId];
|
|
|
|
// detectedPhotons->Fill();
|
|
}
|
|
else
|
|
{
|
|
G4cout << G4endl << "Pointer to detectedPhotons tree is NULL!" << G4endl;
|
|
}
|
|
}
|
|
|
|
char* Analysis::FileName()
|
|
{
|
|
return fileName;
|
|
}
|
|
|
|
void Analysis::SetGpsPosition(G4ThreeVector position)
|
|
{
|
|
gpsPositionX = position[0];
|
|
gpsPositionY = position[1];
|
|
gpsPositionZ = position[2];
|
|
}
|
|
|
|
void Analysis::IncreaseReflectionsAtMirror(Int_t trackId)
|
|
{
|
|
reflectionsAtMirror[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::SetGpsDirection(G4ThreeVector direction)
|
|
{
|
|
gpsDirectionX = direction[0];
|
|
gpsDirectionY = direction[1];
|
|
gpsDirectionZ = direction[2];
|
|
}
|
|
|
|
void Analysis::IncreaseReflectionsAtFibreSurface(Int_t trackId)
|
|
{
|
|
reflectionsAtFibreSurface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseTotalReflectionsAtCladCladInterface(Int_t trackId)
|
|
{
|
|
reflectionsTotalAtCladCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseTotalReflectionsAtCoreCladInterface(Int_t trackId)
|
|
{
|
|
reflectionsTotalAtCoreCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseFresnelReflectionsAtCladCladInterface(Int_t trackId)
|
|
{
|
|
reflectionsFresnelAtCladCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseFresnelReflectionsAtCoreCladInterface(Int_t trackId)
|
|
{
|
|
reflectionsFresnelAtCoreCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseRefractionsAtCladCladInterface(Int_t trackId)
|
|
{
|
|
refractionsAtCladCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseRefractionsAtCoreCladInterface(Int_t trackId)
|
|
{
|
|
refractionsAtCoreCladInterface[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseRayleighScatterings(Int_t trackId)
|
|
{
|
|
rayleighScatterings[trackId]++;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseReflectionRefractionAndScatteringVectors(Int_t trackId)
|
|
{
|
|
for(int i=reflectionsAtMirror.size(); i<trackId+1; i++)
|
|
reflectionsAtMirror.push_back(0);
|
|
|
|
for(int i=reflectionsAtFibreSurface.size(); i<trackId+1; i++)
|
|
reflectionsAtFibreSurface.push_back(0);
|
|
|
|
for(int i=reflectionsTotalAtCladCladInterface.size(); i<trackId+1; i++)
|
|
reflectionsTotalAtCladCladInterface.push_back(0);
|
|
|
|
for(int i=reflectionsTotalAtCoreCladInterface.size(); i<trackId+1; i++)
|
|
reflectionsTotalAtCoreCladInterface.push_back(0);
|
|
|
|
for(int i=reflectionsFresnelAtCladCladInterface.size(); i<trackId+1; i++)
|
|
reflectionsFresnelAtCladCladInterface.push_back(0);
|
|
|
|
for(int i=reflectionsFresnelAtCoreCladInterface.size(); i<trackId+1; i++)
|
|
reflectionsFresnelAtCoreCladInterface.push_back(0);
|
|
|
|
for(int i=refractionsAtCladCladInterface.size(); i<trackId+1; i++)
|
|
refractionsAtCladCladInterface.push_back(0);
|
|
|
|
for(int i=refractionsAtCoreCladInterface.size(); i<trackId+1; i++)
|
|
refractionsAtCoreCladInterface.push_back(0);
|
|
|
|
for(int i=rayleighScatterings.size(); i<trackId+1; i++)
|
|
rayleighScatterings.push_back(0);
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseLengthInCore(Int_t trackId, Float_t lengthValue)
|
|
{
|
|
lengthInCore[trackId] += lengthValue;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseLengthInInnerCladding(Int_t trackId, Float_t lengthValue)
|
|
{
|
|
lengthInInnerCladding[trackId] += lengthValue;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseLengthInOuterCladding(Int_t trackId, Float_t lengthValue)
|
|
{
|
|
lengthInOuterCladding[trackId] += lengthValue;
|
|
}
|
|
|
|
|
|
void Analysis::IncreaseLengthVectors(Int_t trackId)
|
|
{
|
|
for(int i=lengthInCore.size(); i<trackId+1; i++)
|
|
lengthInCore.push_back(0);
|
|
|
|
for(int i=lengthInInnerCladding.size(); i<trackId+1; i++)
|
|
lengthInInnerCladding.push_back(0);
|
|
|
|
for(int i=lengthInOuterCladding.size(); i<trackId+1; i++)
|
|
lengthInOuterCladding.push_back(0);
|
|
}
|