CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Stub.cc
Go to the documentation of this file.
4 
8 
9 #include <cmath>
10 #include <bitset>
11 
12 using namespace std;
13 using namespace trklet;
14 
15 Stub::Stub(Settings const& settings) : settings_(settings) {}
16 
17 Stub::Stub(L1TStub& stub, Settings const& settings, Globals& globals) : settings_(settings) {
18  const string& stubwordhex = stub.stubword();
19 
20  const string stubwordbin = convertHexToBin(stubwordhex);
21 
22  layerdisk_ = stub.layerdisk();
23 
24  int nbendbits = stub.isPSmodule() ? N_BENDBITS_PS : N_BENDBITS_2S;
25 
26  int nalphabits = 0;
27 
28  int nrbits = settings_.nrbitsstub(layerdisk_);
29  int nzbits = settings_.nzbitsstub(layerdisk_);
30  int nphibits = settings_.nphibitsstub(layerdisk_);
31 
32  if (layerdisk_ >= N_LAYER && !stub.isPSmodule()) {
33  nalphabits = settings.nbitsalpha();
34  nrbits = 7;
35  }
36 
37  assert(nbendbits + nalphabits + nrbits + nzbits + nphibits == 36);
38 
39  bitset<32> rbits(stubwordbin.substr(0, nrbits));
40  bitset<32> zbits(stubwordbin.substr(nrbits, nzbits));
41  bitset<32> phibits(stubwordbin.substr(nrbits + nzbits, nphibits));
42  bitset<32> alphabits(stubwordbin.substr(nphibits + nzbits + nrbits, nalphabits));
43  bitset<32> bendbits(stubwordbin.substr(nphibits + nzbits + nrbits + nalphabits, nbendbits));
44 
45  int newbend = bendbits.to_ulong();
46 
47  int newr = rbits.to_ulong();
48  if (layerdisk_ < N_LAYER) {
49  if (newr >= (1 << (nrbits - 1)))
50  newr = newr - (1 << nrbits);
51  }
52 
53  int newz = zbits.to_ulong();
54  if (newz >= (1 << (nzbits - 1)))
55  newz = newz - (1 << nzbits);
56 
57  int newphi = phibits.to_ulong();
58 
59  int newalpha = alphabits.to_ulong();
60  if (newalpha >= (1 << (nalphabits - 1)))
61  newalpha = newalpha - (1 << nalphabits);
62 
63  l1tstub_ = &stub;
64 
65  bend_.set(newbend, nbendbits, true, __LINE__, __FILE__);
66 
67  phi_.set(newphi, nphibits, true, __LINE__, __FILE__);
68  phicorr_.set(newphi, nphibits, true, __LINE__, __FILE__);
69  bool pos = false;
70  if (layerdisk_ >= N_LAYER) {
71  pos = true;
72  int disk = layerdisk_ - N_LAYER + 1;
73  if (stub.z() < 0.0)
74  disk = -disk;
75  disk_.set(disk, 4, false, __LINE__, __FILE__);
76  if (!stub.isPSmodule()) {
77  alpha_.set(newalpha, nalphabits, false, __LINE__, __FILE__);
78  nrbits = 4;
79  }
80  } else {
81  disk_.set(0, 4, false, __LINE__, __FILE__);
82  layer_.set(layerdisk_, 3, true, __LINE__, __FILE__);
83  }
84  r_.set(newr, nrbits, pos, __LINE__, __FILE__);
85  z_.set(newz, nzbits, false, __LINE__, __FILE__);
86 
87  if (settings.writeMonitorData("StubBend")) {
88  unsigned int nsimtrks = globals.event()->nsimtracks();
89 
90  for (unsigned int isimtrk = 0; isimtrk < nsimtrks; isimtrk++) {
91  const L1SimTrack& simtrk = globals.event()->simtrack(isimtrk);
92  if (stub.tpmatch2(simtrk.trackid())) {
93  double dr = 0.18;
94  double rinv = simtrk.charge() * 0.01 * settings_.c() * settings_.bfield() / simtrk.pt();
95  double pitch = settings_.stripPitch(stub.isPSmodule());
96  double bend = stub.r() * dr * 0.5 * rinv / pitch;
97 
98  globals.ofstream("stubbend.dat") << layerdisk_ << " " << stub.isPSmodule() << " "
99  << simtrk.pt() * simtrk.charge() << " " << bend << " " << newbend << " "
100  << settings.benddecode(newbend, layerdisk_, stub.isPSmodule()) << " "
101  << settings.bendcut(newbend, layerdisk_, stub.isPSmodule()) << endl;
102  }
103  }
104  }
105 }
106 
107 FPGAWord Stub::iphivmFineBins(int VMbits, int finebits) const {
108  unsigned int finephi = (phicorr_.value() >> (phicorr_.nbits() - VMbits - finebits)) & ((1 << finebits) - 1);
109  return FPGAWord(finephi, finebits, true, __LINE__, __FILE__);
110 }
111 
112 unsigned int Stub::phiregionaddress() const {
114  return (iphi << 7) + stubindex_.value();
115 }
116 
119  FPGAWord phiregion(iphi, 3, true, __LINE__, __FILE__);
120  return phiregion.str() + stubindex_.str();
121 }
122 
123 void Stub::setAllStubIndex(int nstub) {
124  if (nstub >= (1 << N_BITSMEMADDRESS)) {
125  if (settings_.debugTracklet())
126  edm::LogPrint("Tracklet") << "Warning too large stubindex!";
127  nstub = (1 << N_BITSMEMADDRESS) - 1;
128  }
129 
131 }
132 
133 void Stub::setPhiCorr(int phiCorr) {
134  int iphicorr = phi_.value() - phiCorr;
135 
136  if (iphicorr < 0)
137  iphicorr = 0;
138  if (iphicorr >= (1 << phi_.nbits()))
139  iphicorr = (1 << phi_.nbits()) - 1;
140 
141  phicorr_.set(iphicorr, phi_.nbits(), true, __LINE__, __FILE__);
142 }
143 
144 double Stub::rapprox() const {
145  if (disk_.value() == 0) {
146  int lr = 1 << (8 - settings_.nrbitsstub(layer_.value()));
147  return r_.value() * settings_.kr() * lr + settings_.rmean(layer_.value());
148  }
149  if (!l1tstub_->isPSmodule()) {
150  if (abs(disk_.value()) <= 2)
151  return settings_.rDSSinner(r_.value());
152  else
153  return settings_.rDSSouter(r_.value());
154  }
155  return r_.value() * settings_.kr();
156 }
157 
158 double Stub::zapprox() const {
159  if (disk_.value() == 0) {
160  int lz = 1;
161  if (layer_.value() >= 3) {
162  lz = 16;
163  }
164  return z_.value() * settings_.kz() * lz;
165  }
166  int sign = 1;
167  if (disk_.value() < 0)
168  sign = -1;
169  if (sign < 0) {
170  //Should understand why this is needed to get agreement with integer calculations
171  return (z_.value() + 1) * settings_.kz() + sign * settings_.zmean(abs(disk_.value()) - 1);
172  } else {
173  return z_.value() * settings_.kz() + sign * settings_.zmean(abs(disk_.value()) - 1);
174  }
175 }
176 
177 double Stub::phiapprox(double phimin, double) const {
178  int lphi = 1;
179  if (layer_.value() >= 3) {
180  lphi = 8;
181  }
182  return reco::reduceRange(phimin + phi_.value() * settings_.kphi() / lphi);
183 }
184 
185 unsigned int Stub::layerdisk() const {
186  if (layer_.value() == -1)
187  return N_LAYER - 1 + abs(disk_.value());
188  return layer_.value();
189 }
trklet::Stub::phi_
FPGAWord phi_
Definition: Stub.h:90
trklet::Stub::Stub
Stub(Settings const &settings)
Definition: Stub.cc:15
MessageLogger.h
trklet::Stub::phiregionaddressstr
std::string phiregionaddressstr() const
Definition: Stub.cc:117
trklet::FPGAWord::str
std::string str() const
Definition: FPGAWord.cc:54
deltaPhi.h
phimin
float phimin
Definition: ReggeGribovPartonMCHadronizer.h:107
trklet::Stub::bend
const FPGAWord & bend() const
Definition: Stub.h:58
trklet::Stub::r_
FPGAWord r_
Definition: Stub.h:88
trklet::Stub::layer_
FPGAWord layer_
Definition: Stub.h:86
pos
Definition: PixelAliasList.h:18
trklet::Settings
Definition: Settings.h:52
trklet::L1TStub
Definition: L1TStub.h:14
trklet::L1TStub::z
double z() const
Definition: L1TStub.h:57
trklet::convertHexToBin
std::string convertHexToBin(const std::string &stubwordhex)
Definition: Util.h:60
trklet::Settings::nbitsallstubs
unsigned int nbitsallstubs(unsigned int layerdisk) const
Definition: Settings.h:106
cms::cuda::assert
assert(be >=bs)
trklet::FPGAWord::set
void set(int value, int nbits, bool positive=true, int line=-1, const char *file=nullptr)
Definition: FPGAWord.cc:14
Validation_hcalonly_cfi.sign
sign
Definition: Validation_hcalonly_cfi.py:32
trklet::L1SimTrack
Definition: L1SimTrack.h:14
trklet::Settings::rmean
double rmean(unsigned int iLayer) const
Definition: Settings.h:164
trklet::Stub::phiapprox
double phiapprox(double phimin, double) const
Definition: Stub.cc:177
trklet::FPGAWord::nbits
int nbits() const
Definition: FPGAWord.h:25
trklet::Stub::bend_
FPGAWord bend_
Definition: Stub.h:93
trklet::Settings::rDSSinner
double rDSSinner(unsigned int iBin) const
Definition: Settings.h:171
trklet::Globals
Definition: Globals.h:30
trklet::Stub::phicorr_
FPGAWord phicorr_
Definition: Stub.h:95
trklet::Settings::bfield
double bfield() const
Definition: Settings.h:253
trklet::Settings::nrbitsstub
unsigned int nrbitsstub(unsigned int layerdisk) const
Definition: Settings.h:84
trklet::Stub::disk_
FPGAWord disk_
Definition: Stub.h:87
trklet::Settings::rDSSouter
double rDSSouter(unsigned int iBin) const
Definition: Settings.h:174
LEDCalibrationChannels.iphi
iphi
Definition: LEDCalibrationChannels.py:64
trklet::Stub::z_
FPGAWord z_
Definition: Stub.h:89
trklet::L1SimTrack::trackid
int trackid() const
Definition: L1SimTrack.h:24
trklet::Settings::benddecode
double benddecode(int ibend, int layerdisk, bool isPSmodule) const
Definition: Settings.h:402
trklet::Stub::disk
const FPGAWord & disk() const
Definition: Stub.h:68
trklet::Stub::iphivmFineBins
FPGAWord iphivmFineBins(int VMbits, int finebits) const
Definition: Stub.cc:107
trklet::N_BENDBITS_2S
constexpr unsigned int N_BENDBITS_2S
Definition: Settings.h:30
trklet::L1TStub::isPSmodule
unsigned int isPSmodule() const
Definition: L1TStub.h:94
trklet::Settings::bendcut
double bendcut(int ibend, int layerdisk, bool isPSmodule) const
Definition: Settings.h:410
trklet::N_LAYER
constexpr int N_LAYER
Definition: Settings.h:21
trklet::Settings::zmean
double zmean(unsigned int iDisk) const
Definition: Settings.h:167
SLHCEvent.h
trklet::Stub::layerdisk
unsigned int layerdisk() const
Definition: Stub.cc:185
trklet::L1TStub::layerdisk
int layerdisk() const
Definition: L1TStub.h:104
trklet::FPGAWord
Definition: FPGAWord.h:9
trklet::Settings::kr
double kr() const
Definition: Settings.h:304
trklet::L1TStub::tpmatch2
bool tpmatch2(int tp) const
Definition: L1TStub.cc:111
trklet::rinv
double rinv(double phi1, double phi2, double r1, double r2)
Definition: Util.h:49
trklet::L1SimTrack::pt
double pt() const
Definition: L1SimTrack.h:26
trklet::Settings::nzbitsstub
unsigned int nzbitsstub(unsigned int layerdisk) const
Definition: Settings.h:82
trklet::Stub::rapprox
double rapprox() const
Definition: Stub.cc:144
trklet::SLHCEvent::nsimtracks
unsigned int nsimtracks() const
Definition: SLHCEvent.h:59
Globals.h
trklet::Settings::nphibitsstub
unsigned int nphibitsstub(unsigned int layerdisk) const
Definition: Settings.h:83
trklet
Definition: AllInnerStubsMemory.h:10
trklet::Stub::setAllStubIndex
void setAllStubIndex(int nstub)
Definition: Stub.cc:123
trklet::FPGAWord::value
int value() const
Definition: FPGAWord.h:24
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
trklet::Settings::nbitsalpha
int nbitsalpha() const
Definition: Settings.h:217
trklet::Settings::stripPitch
double stripPitch(bool isPSmodule) const
Definition: Settings.h:260
trklet::Settings::kphi
double kphi() const
Definition: Settings.h:298
trklet::L1SimTrack::charge
int charge() const
Definition: L1SimTrack.h:34
trklet::N_BITSMEMADDRESS
constexpr unsigned int N_BITSMEMADDRESS
Definition: Settings.h:39
trklet::Settings::writeMonitorData
bool writeMonitorData(std::string module) const
Definition: Settings.h:109
trklet::Stub::phiregionaddress
unsigned int phiregionaddress() const
Definition: Stub.cc:112
std
Definition: JetResolutionObject.h:76
trklet::Stub::setPhiCorr
void setPhiCorr(int phiCorr)
Definition: Stub.cc:133
trklet::L1TStub::stubword
const std::string & stubword() const
Definition: L1TStub.h:108
trklet::Globals::ofstream
std::ofstream & ofstream(std::string fname)
Definition: Globals.cc:44
trklet::Globals::event
SLHCEvent *& event()
Definition: Globals.h:36
flavorHistoryFilter_cfi.dr
dr
Definition: flavorHistoryFilter_cfi.py:37
trklet::Settings::kz
double kz() const
Definition: Settings.h:302
trklet::Settings::debugTracklet
bool debugTracklet() const
Definition: Settings.h:182
Exception.h
trklet::L1TStub::r
double r() const
Definition: L1TStub.h:58
trklet::N_BENDBITS_PS
constexpr unsigned int N_BENDBITS_PS
Definition: Settings.h:29
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
trklet::Stub::layerdisk_
unsigned int layerdisk_
Definition: Stub.h:84
trklet::SLHCEvent::simtrack
const L1SimTrack & simtrack(int i) const
Definition: SLHCEvent.h:61
trklet::Stub::settings_
Settings const & settings_
Definition: Stub.h:100
trklet::Settings::c
double c() const
Definition: Settings.h:212
edm::Log
Definition: MessageLogger.h:70
Stub.h
trklet::Stub::stubindex_
FPGAWord stubindex_
Definition: Stub.h:97
trklet::Stub::zapprox
double zapprox() const
Definition: Stub.cc:158
trklet::Stub::l1tstub_
L1TStub * l1tstub_
Definition: Stub.h:99
trklet::Stub::alpha_
FPGAWord alpha_
Definition: Stub.h:91
reco::reduceRange
constexpr T reduceRange(T x)
Definition: deltaPhi.h:18