CMS 3D CMS Logo

BeamSpotOnlineProducer.cc
Go to the documentation of this file.
1 
4 
10 
12 
13 
14 using namespace edm;
15 
16 
18  : changeFrame_ ( iconf.getParameter<bool> ("changeToCMSCoordinates") )
19  , theMaxZ ( iconf.getParameter<double>("maxZ") )
20  , theSetSigmaZ ( iconf.getParameter<double>("setSigmaZ") )
21  , scalerToken_ ( consumes<BeamSpotOnlineCollection> ( iconf.getParameter<InputTag>("src") ) )
22  , l1GtEvmReadoutRecordToken_ ( consumes<L1GlobalTriggerEvmReadoutRecord>( iconf.getParameter<InputTag>("gtEvmLabel")) )
23  , theBeamShoutMode ( iconf.getUntrackedParameter<unsigned int> ("beamMode",11) )
24 {
25 
26  theMaxR2 = iconf.getParameter<double>("maxRadius");
28 
29  produces<reco::BeamSpot>();
30 
31 }
32 
34 
35 void
37 {
38  //shout MODE only in stable beam
39  bool shoutMODE=false;
41  if (iEvent.getByToken(l1GtEvmReadoutRecordToken_, gtEvmReadoutRecord)){
42  if (gtEvmReadoutRecord->gtfeWord().beamMode() == theBeamShoutMode) shoutMODE=true;
43  }
44  else{
45  shoutMODE=true;
46  }
47 
48  // get scalar collection
50  iEvent.getByToken( scalerToken_, handleScaler);
51 
52  // beam spot scalar object
53  BeamSpotOnline spotOnline;
54 
55  // product is a reco::BeamSpot object
56  auto result = std::make_unique<reco::BeamSpot>();
57 
58  reco::BeamSpot aSpot;
59 
60  bool fallBackToDB=false;
61  if (!handleScaler->empty()){
62  // get one element
63  spotOnline = * ( handleScaler->begin() );
64 
65  // in case we need to switch to LHC reference frame
66  // ignore for the moment rotations, and translations
67  double f = 1.;
68  if (changeFrame_) f = -1.;
69 
70  reco::BeamSpot::Point apoint( f* spotOnline.x(), spotOnline.y(), f* spotOnline.z() );
71 
73  matrix(0,0) = spotOnline.err_x()*spotOnline.err_x();
74  matrix(1,1) = spotOnline.err_y()*spotOnline.err_y();
75  matrix(2,2) = spotOnline.err_z()*spotOnline.err_z();
76  matrix(3,3) = spotOnline.err_sigma_z()*spotOnline.err_sigma_z();
77 
78  double sigmaZ = spotOnline.sigma_z();
79  if (theSetSigmaZ>0)
80  sigmaZ = theSetSigmaZ;
81 
82  aSpot = reco::BeamSpot( apoint,
83  sigmaZ,
84  spotOnline.dxdz(),
85  f* spotOnline.dydz(),
86  spotOnline.width_x(),
87  matrix);
88 
89  aSpot.setBeamWidthY( spotOnline.width_y() );
90  aSpot.setEmittanceX( 0. );
91  aSpot.setEmittanceY( 0. );
92  aSpot.setbetaStar( 0.) ;
93  aSpot.setType( reco::BeamSpot::LHC ); // flag value from scalars
94 
95  // check if we have a valid beam spot fit result from online DQM
96  if ( spotOnline.x() == 0 &&
97  spotOnline.y() == 0 &&
98  spotOnline.z() == 0 &&
99  spotOnline.width_x() == 0 &&
100  spotOnline.width_y() == 0 )
101  {
102  if (shoutMODE){
103  edm::LogWarning("BeamSpotFromDB")
104  << "Online Beam Spot producer falls back to DB value because the scaler values are zero ";
105  }
106  fallBackToDB=true;
107  }
108  double r2=spotOnline.x()*spotOnline.x() + spotOnline.y()*spotOnline.y();
109  if (fabs(spotOnline.z())>=theMaxZ || r2>=theMaxR2){
110  if (shoutMODE){
111  edm::LogError("BeamSpotFromDB")
112  << "Online Beam Spot producer falls back to DB value because the scaler values are too big to be true :"
113  <<spotOnline.x()<<" "<<spotOnline.y()<<" "<<spotOnline.z();
114  }
115  fallBackToDB=true;
116  }
117  }
118  else{
119  //empty online beamspot collection: FED data was empty
120  //the error should probably have been send at unpacker level
121  fallBackToDB=true;
122  }
123 
124  if (fallBackToDB){
125 
127  iSetup.get<BeamSpotObjectsRcd>().get(beamhandle);
128  const BeamSpotObjects *spotDB = beamhandle.product();
129 
130  // translate from BeamSpotObjects to reco::BeamSpot
131  reco::BeamSpot::Point apoint( spotDB->GetX(), spotDB->GetY(), spotDB->GetZ() );
132 
134  for ( int i=0; i<7; ++i ) {
135  for ( int j=0; j<7; ++j ) {
136  matrix(i,j) = spotDB->GetCovariance(i,j);
137  }
138  }
139 
140  // this assume beam width same in x and y
141  aSpot = reco::BeamSpot( apoint,
142  spotDB->GetSigmaZ(),
143  spotDB->Getdxdz(),
144  spotDB->Getdydz(),
145  spotDB->GetBeamWidthX(),
146  matrix );
147  aSpot.setBeamWidthY( spotDB->GetBeamWidthY() );
148  aSpot.setEmittanceX( spotDB->GetEmittanceX() );
149  aSpot.setEmittanceY( spotDB->GetEmittanceY() );
150  aSpot.setbetaStar( spotDB->GetBetaStar() );
151  aSpot.setType( reco::BeamSpot::Tracker );
152 
153  }
154 
155  *result = aSpot;
156 
157  iEvent.put(std::move(result));
158 
159 }
160 
double Getdydz() const
get dydz slope, crossing angle in YZ
std::vector< BeamSpotOnline > BeamSpotOnlineCollection
math::Error< dimension >::type CovarianceMatrix
Definition: BeamSpot.h:31
T getParameter(std::string const &) const
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:136
float err_y() const
float dxdz() const
double GetY() const
get Y beam position
const edm::EDGetTokenT< L1GlobalTriggerEvmReadoutRecord > l1GtEvmReadoutRecordToken_
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:519
float err_x() const
float dydz() const
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
void produce(edm::Event &iEvent, const edm::EventSetup &iSetup) override
produce a beam spot class
const L1GtfeExtWord gtfeWord() const
get / set GTFE word (record) in the GT readout record
BeamSpotOnlineProducer(const edm::ParameterSet &iConf)
constructor
float sigma_z() const
const edm::EDGetTokenT< BeamSpotOnlineCollection > scalerToken_
float y() const
float x() const
double GetSigmaZ() const
get sigma Z, RMS bunch length
double GetBeamWidthX() const
get average transverse beam width
math::XYZPoint Point
point in the space
Definition: BeamSpot.h:29
float z() const
float width_y() const
double GetBeamWidthY() const
get average transverse beam width
double GetEmittanceX() const
get emittance
int iEvent
Definition: GenABIO.cc:230
double GetZ() const
get Z beam position
double f[11][100]
const unsigned int theBeamShoutMode
~BeamSpotOnlineProducer() override
destructor
double Getdxdz() const
get dxdz slope, crossing angle in XZ
float err_z() const
double GetX() const
get X beam position
const T & get() const
Definition: EventSetup.h:59
double GetBetaStar() const
get beta star
double GetCovariance(int i, int j) const
get i,j element of the full covariance matrix 7x7
HLT enums.
float width_x() const
double GetEmittanceY() const
get emittance
float err_sigma_z() const
const cms_uint16_t beamMode() const
T const * product() const
Definition: ESHandle.h:86
def move(src, dest)
Definition: eostools.py:510