Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:19:00

0001 #ifndef IOMC_HLLHCEvtVtxGenerator_H
0002 #define IOMC_HLLHCEvtVtxGenerator_H
0003 
0004 /**
0005  * Generate event vertices given beams sizes, crossing angle
0006  * offset, and crab rotation. 
0007  * Attention: All values are assumed to be mm for spatial coordinates
0008  * and ns for time.
0009  * Attention: This class fix the the vertex time generation of HLLHCEvtVtxGenerator
0010  *
0011  * $Id: HLLHCEvtVtxGenerator_Fix.h,v 1.0 2015/03/15 10:34:38 Exp $
0012  */
0013 
0014 #include "IOMC/EventVertexGenerators/interface/BaseEvtVtxGenerator.h"
0015 
0016 #include <string>
0017 
0018 namespace CLHEP {
0019   class RandFlat;
0020 }
0021 
0022 namespace edm {
0023   class ConfigurationDescriptions;
0024 }
0025 
0026 class HLLHCEvtVtxGenerator : public BaseEvtVtxGenerator {
0027 public:
0028   HLLHCEvtVtxGenerator(const edm::ParameterSet& p);
0029 
0030   /** Copy constructor */
0031   HLLHCEvtVtxGenerator(const HLLHCEvtVtxGenerator& p) = delete;
0032 
0033   /** Copy assignment operator */
0034   HLLHCEvtVtxGenerator& operator=(const HLLHCEvtVtxGenerator& rhs) = delete;
0035 
0036   ~HLLHCEvtVtxGenerator() override;
0037 
0038   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0039 
0040   /// return a new event vertex
0041   HepMC::FourVector newVertex(CLHEP::HepRandomEngine*) const override;
0042 
0043   TMatrixD const* GetInvLorentzBoost() const override { return nullptr; };
0044 
0045 private:
0046   //spatial and time offset for mean collision
0047   const double fMeanX, fMeanY, fMeanZ, fTimeOffset;
0048 
0049   //proton beam energy
0050   const double momeV;
0051   const double gamma;
0052   const double beta;
0053   const double betagamma;
0054 
0055   //crossing angle
0056   const double phi;
0057 
0058   //crab cavity frequency
0059   const double wcc;
0060 
0061   // 800 MHz RF?
0062   const bool RF800;
0063 
0064   //beta crossing plane (m)
0065   const double betx;
0066 
0067   //beta separation plane (m)
0068   const double bets;
0069 
0070   //horizontal emittance
0071   const double epsxn;
0072 
0073   //vertical emittance
0074   const double epssn;
0075 
0076   //bunch length
0077   const double sigs;
0078 
0079   //crabbing angle crossing
0080   const double alphax;
0081 
0082   //crabbing angle separation
0083   const double alphay;
0084 
0085   // ratio of crabbing angle to crossing angle
0086   const double oncc;
0087 
0088   //normalized crossing emittance
0089   const double epsx;
0090 
0091   //normlaized separation emittance
0092   const double epss;
0093 
0094   //size in x
0095   const double sigx;
0096 
0097   // crossing angle * crab frequency
0098   const double phiCR;
0099 
0100   //width for y plane
0101   double sigma(double z, double epsilon, double beta, double betagamma) const;
0102 
0103   //density with crabbing
0104   double integrandCC(double x, double z, double t) const;
0105 
0106   // 4D intensity
0107   double intensity(double x, double y, double z, double t) const;
0108 };
0109 
0110 #endif