//////////////////////////////////////////////////////////////////// /// \class RAT::VertexGen_ELLIE /// /// \brief Vertex generator - for ELLIE (or other light source) beams /// /// \author Jeanne Wilson -- contact person /// /// REVISION HISTORY /// - first version : J. Wilson - a modification of Gen_ELLIE to only deal with vertex. /// motivated by need for separate vertex, position generators for /// use with coincidence generator /// - 2014-11-26: Matt Strait - doxygen fixes /// - 2017-02-07: Nuno Barros - Added BeginOfRun to load the database /// - 2021-02-25: Jeanne Wilson, Esther Turner, Daniel Cookman - Adding 2D beam angular distributions /// /// \details This generator simulates bunches of photons from fibres. User RatDB settings for fibre position, direction, angular, timing and wavelength distributions. //////////////////////////////////////////////////////////////////// #ifndef __RAT_VertexGen_ELLIE__ #define __RAT_VertexGen_ELLIE__ #include #include #include #include #include #include "RAT/GLG4Gen.hh" #include class G4ParticleDefinition; namespace RAT { class VertexGen_ELLIE : public GLG4VertexGen { public: VertexGen_ELLIE(const char *arg_dbname="vertexellie"); virtual ~VertexGen_ELLIE(); void Init(); virtual void BeginOfRun(); virtual void GeneratePrimaryVertex(G4Event *event, G4ThreeVector &dx, G4double dt); float SampleRandom(std::vector myx, std::vector myy, std::vector mycum); /** interface for getting and setting the state*/ virtual void SetState(G4String state); virtual G4String GetState(); private: //---------Initialisation methods---------------------- /** Initialise intensity info from ELLIE DB, noting photon thinning & whether to Poisson fluctuate. */ void InitialiseIntensity(DBLinkPtr lELLIE); /** Initialise ELLIE fibre info from DB, returning a link to the fibre DB entry.*/ DBLinkPtr InitialiseFibreInfo(DBLinkPtr lELLIE); /** Initialise source info from DB, returning a link to the source info DB entry. */ DBLinkPtr InitialiseSourceInfo(DBLinkPtr lELLIE); /** Initialise wavelength info from DB. */ void InitialiseWavelengthInfo(DBLinkPtr lELLIE, DBLinkPtr lsource); /** Intialise wavelength distribution info from DB. */ void InitialiseWavelengthDistribution(DBLinkPtr lELLIE, DBLinkPtr lsource); /** Initialise wavelength distribution string, returning a link to the distribution DB info entry.*/ DBLinkPtr InitialiseWavelengthDistString(DBLinkPtr lELLIE, DBLinkPtr lsource); /** Initialise wavelength distribution x-values from DB. */ void InitialiseWavelengthDistXVals(DBLinkPtr lELLIEwave); /** Initialise wavelength distribution y-values from DB. */ void InitialiseWavelengthDistYVals(DBLinkPtr lELLIEwave); /** Initialise a subset of the wavelength distribution if a range is specified. */ void InitWavelengthDistSubset(DBLinkPtr lELLIE); /** Initialise timing info from DB. */ void InitialiseTimingInfo(DBLinkPtr lELLIE, DBLinkPtr lsource); /** Initialise angular info from DB. */ void InitialiseAngularInfo(DBLinkPtr lELLIE, DBLinkPtr lfibre); /** Initialise 1D angular distribution. */ void Initialise1DAngularDist(); /** Initialise 2D angular distribution. */ void Initialise2DAngularDist(DBLinkPtr lELLIE); //---------Vertex generation methods------------------ void FindThreeNearestNeighbours(double theta, double phi); G4ThreeVector Find_Projection(G4ThreeVector& vectorToBeProjected, G4ThreeVector& normalToPlane); //Finding the distance between the two input points double FindDistance(G4ThreeVector& photonHitPoint, G4ThreeVector& thisPMT); //Doing an an inverse distance weighted average of the three nearest neighbours double Interpolation(); void SampleRandom2D(double& theta, double& phi); //------Private members----------------------------------- /** Only call the initialisation once - if this flag isn't set */ bool fInitDB; /** Fibre IDs */ std::string fFibreID; /** Source */ std::string fSource; /** ELLIE direction */ double fDirX, fDirY, fDirZ; /** ELLIE position */ double fPosX, fPosY, fPosZ; /** Wavelength mode for this ELLIE (1 for "mono", 0 for the distribution) The distribution name is the wavelength as a string and defines the index of the ELLIEWAVE ratdb table*/ int fMonoWls; /** ELLIE wavelength distribution - string tells you what distribution to read */ std::string fELLIEWavelength; /** ELLIE wavelength range - use a portion of the wavelength distribution within low and high setpoints */ std::vector fELLIEWavelengthRange; /** Time mode for this ELLIE (1 for "mono", 0 for the distribution) The distribution name is the wavelength as a string and defines the index of the ELLIETIME ratdb table*/ int fMonoTime; /** ELLIE time distribution - string tells you what distribution to read */ std::string fELLIETime; /** Angle mode for this ELLIE (1 for "iso", 0 otherwise) The distribution name is the fibre id and defines the index of the ELLIEANGLE ratdb table*/ int fIsoAngle; /** flag set to 1 if there is to be no dispersion of the beam */ int fZero_angle; /** ELLIE angular distribution - string tells you what distribution to read */ std::string fELLIEAngle; /** Boolean flag that determines if a 1D or 2D distribution is being used. false for a 1D distribution, true for a 2D distribution.*/ bool f2DAngle; /** true = number of photons per pulse is selected from poisson distribution with mean=photons_per_event, \n false = number of photons per pulse=photons_per_event */ bool fPoisson; /** Vector for specified angles */ std::vector fAng; /** Vector for cumulative angular intensity distribution */ std::vector fAngCumMag; /** Vector for angular intensity */ std::vector fAngInt; /** Or for 2dimensional angular distribution */ std::vector fAngPmts; // list of PMTs within angular cut of the beamspot centre std::vector fAngProbs; // relative intensity at each of the PMTs //std::vector fAng_errors; /** Vectors for PMT info */ std::vector fAllPMTXs; std::vector fAllPMTYs; std::vector fAllPMTZs; /** Vector for the position for just the PMTs that fall within our angle cut */ std::vector fPmtPositions; /** Vectors for thetas and phis of these PMTs in the fibre coordinates */ std::vector fPmtPhis; std::vector fPmtThetas; /** Vector for specified times */ std::vector fTime; /** Vector for cumulative time intensity distribution */ std::vector fTimeCumMag; /** Vector for time intensity */ std::vector fTimeInt; /** Vector for specified wavelengths */ std::vector fWave; /** Vector for cumulative wave intensity distribution */ std::vector fWaveCumMag; /** Vector for wave intensity */ std::vector fWaveInt; /** Factor for converting angles (random generator uses radians so if distribution set in degrees = 360/twopi, otherwise =1*/ float fConvertAngle; /** Number of photons per ELLIE pulse */ int fPhotonsPerEvent; /** Particles must be defined as optical photons */ G4ParticleDefinition *fPhotonDef; //** Direction, position and perpendicular to direction vectors for the fibre */ G4ThreeVector fDir; G4ThreeVector fPerp; G4ThreeVector fPos; G4ThreeVector fNormal; G4ThreeVector fZDirProj; //Position of the centre of beam on the opposite side of the PSUP G4ThreeVector fBeamCentre; //** Variables for the importance sampling */ double fMaxTheta; //The maximum theta of any of the PMTs in the profile double fMean; //The mean of the Gaussian used as the envelope double fStdDev; //The standard deviation of the Gaussian used as the envelope double fAmp; //The amplitude of the Gaussian used as the envelope //Vectors of the distances to the three nearest PMTs to an abitrary point and the intensity observed in those PMTs std::vector fDistances; std::vector fIntensities; }; } // namespace RAT #endif