File indexing completed on 2024-04-06 12:15:19
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include <memory>
0021 #include <fstream>
0022 #include <iomanip>
0023
0024
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0027 #include "FWCore/Framework/interface/Event.h"
0028 #include "FWCore/Framework/interface/MakerMacros.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0030
0031 #include "SimDataFormats/TrackingHit/interface/PSimHitContainer.h"
0032 #include "SimDataFormats/TrackingHit/interface/PSimHit.h"
0033 #include <DataFormats/RPCDigi/interface/RPCDigi.h>
0034 #include <DataFormats/RPCDigi/interface/RPCDigiCollection.h>
0035 #include <DataFormats/MuonDetId/interface/RPCDetId.h>
0036
0037 #include <Geometry/Records/interface/MuonGeometryRecord.h>
0038 #include <Geometry/RPCGeometry/interface/RPCGeometry.h>
0039 #include <Geometry/RPCGeometry/interface/RPCGeomServ.h>
0040 #include <Geometry/CommonTopologies/interface/RectangularStripTopology.h>
0041 #include <Geometry/CommonTopologies/interface/TrapezoidalStripTopology.h>
0042
0043 class RPCGEO : public edm::one::EDAnalyzer<> {
0044 public:
0045 explicit RPCGEO(const edm::ParameterSet&);
0046 ~RPCGEO() override;
0047
0048 void beginJob() override {}
0049 void analyze(edm::Event const& iEvent, edm::EventSetup const&) override;
0050 void endJob() override {}
0051
0052 private:
0053 const edm::ESGetToken<RPCGeometry, MuonGeometryRecord> tokRPC_;
0054 std::ofstream detidsstream;
0055 std::ofstream detailstream;
0056 };
0057
0058 RPCGEO::RPCGEO(const edm::ParameterSet& )
0059 : tokRPC_{esConsumes<RPCGeometry, MuonGeometryRecord>(edm::ESInputTag{})} {
0060 detailstream.open("RPCGeometry.out");
0061 detidsstream.open("RPCDetIdLst.out");
0062 edm::LogVerbatim("RPCGeometry") << "Opening output files :: RPCGeometry.out and RPCDetIdLst.out";
0063 }
0064
0065 RPCGEO::~RPCGEO() {
0066 detailstream.close();
0067 detidsstream.close();
0068 edm::LogVerbatim("RPCGeometry") << "Closing output files :: RPCGeometry.out and RPCDetIdLst.out";
0069 }
0070
0071
0072 void RPCGEO::analyze(const edm::Event& , const edm::EventSetup& iSetup) {
0073 using namespace edm;
0074
0075 detailstream << " Getting the RPC Geometry" << std::endl;
0076 const RPCGeometry* rpcGeo = &iSetup.getData(tokRPC_);
0077
0078 int StripsInCMS = 0;
0079 int counterstripsBarrel = 0;
0080 int counterstripsEndCap = 0;
0081 int RollsInCMS = 0;
0082 int counterRollsBarrel = 0;
0083 int counterRollsMB1MB2MB3 = 0;
0084 int counterRollsMB4 = 0;
0085 int RB4Wm2 = 0;
0086 int RB4Wm1 = 0;
0087 int RB4W0 = 0;
0088 int RB4W1 = 0;
0089 int RB4W2 = 0;
0090 int counterRollsEndCap = 0;
0091 int ENDCAP[5][4];
0092 int ENDCAProll[5][4];
0093 int rollsNearDiskp3 = 0;
0094 int rollsNearDiskp2 = 0;
0095 float sumstripwbarrel = 0;
0096 float sumstripwendcap = 0;
0097 float areabarrel = 0;
0098 float areaendcap = 0;
0099
0100
0101
0102
0103 for (int i = 1; i < 5; i++) {
0104 for (int j = 1; j < 4; j++) {
0105 ENDCAP[i][j] = 0;
0106 }
0107 }
0108
0109 for (int i = 1; i < 5; i++) {
0110 for (int j = 1; j < 4; j++) {
0111 ENDCAProll[i][j] = 0;
0112 }
0113 }
0114
0115 for (TrackingGeometry::DetContainer::const_iterator it = rpcGeo->dets().begin(); it < rpcGeo->dets().end(); it++) {
0116
0117
0118
0119 if (dynamic_cast<const RPCChamber*>(*it) != nullptr) {
0120 const RPCChamber* ch = dynamic_cast<const RPCChamber*>(*it);
0121 std::vector<const RPCRoll*> rolls = (ch->rolls());
0122
0123 if (ch->id().region() == 1) {
0124 switch (ch->id().station()) {
0125 case 1:
0126 switch (ch->id().ring()) {
0127 case 1:
0128 ++ENDCAP[ch->id().station()][ch->id().ring()];
0129 break;
0130 case 2:
0131 ++ENDCAP[ch->id().station()][ch->id().ring()];
0132 break;
0133 case 3:
0134 ++ENDCAP[ch->id().station()][ch->id().ring()];
0135 break;
0136 }
0137 break;
0138 case 2:
0139 switch (ch->id().ring()) {
0140 case 1:
0141 ++ENDCAP[ch->id().station()][ch->id().ring()];
0142 break;
0143 case 2:
0144 ++ENDCAP[ch->id().station()][ch->id().ring()];
0145 break;
0146 case 3:
0147 ++ENDCAP[ch->id().station()][ch->id().ring()];
0148 break;
0149 }
0150 break;
0151 case 3:
0152 switch (ch->id().ring()) {
0153 case 1:
0154 ++ENDCAP[ch->id().station()][ch->id().ring()];
0155 break;
0156 case 2:
0157 ++ENDCAP[ch->id().station()][ch->id().ring()];
0158 break;
0159 case 3:
0160 ++ENDCAP[ch->id().station()][ch->id().ring()];
0161 break;
0162 }
0163 break;
0164 case 4:
0165 switch (ch->id().ring()) {
0166 case 1:
0167 ++ENDCAP[ch->id().station()][ch->id().ring()];
0168 break;
0169 case 2:
0170 ++ENDCAP[ch->id().station()][ch->id().ring()];
0171 break;
0172 case 3:
0173 ++ENDCAP[ch->id().station()][ch->id().ring()];
0174 break;
0175 }
0176 break;
0177 }
0178 }
0179
0180
0181
0182 for (std::vector<const RPCRoll*>::const_iterator r = rolls.begin(); r != rolls.end(); ++r) {
0183 RPCDetId rpcId = (*r)->id();
0184 int stripsinthisroll = (*r)->nstrips();
0185 RPCGeomServ rpcsrv(rpcId);
0186 ++RollsInCMS;
0187
0188
0189
0190
0191
0192
0193 if (rpcId.region() == 0) {
0194
0195 const RectangularStripTopology* top_ = dynamic_cast<const RectangularStripTopology*>(&((*r)->topology()));
0196 float s1 = static_cast<float>(1) - 0.5;
0197 float sLast = static_cast<float>(stripsinthisroll) - 0.5;
0198 float stripl = top_->stripLength();
0199 float stripw = top_->pitch();
0200 areabarrel = areabarrel + stripl * stripw * stripsinthisroll;
0201 sumstripwbarrel = sumstripwbarrel + stripw * stripsinthisroll;
0202 std::string name = rpcsrv.name();
0203
0204
0205
0206 detidsstream << rpcId.rawId() << " \t " << std::setw(24) << name << " \t " << rpcId << std::endl;
0207 detailstream << "All Info " << rpcId.rawId() << " = " << std::setw(24) << name << " || ";
0208 detailstream << " strips: length =" << stripl << "[cm] & width =" << stripw
0209 << "[cm] & number =" << stripsinthisroll;
0210 detailstream << " || area roll =" << stripl * stripw << "[cm^2] || area total barrel = " << areabarrel
0211 << "[cm^2]" << std::endl;
0212
0213
0214 const BoundPlane& RPCSurface = (*r)->surface();
0215 GlobalPoint FirstStripCenterPointInGlobal = RPCSurface.toGlobal(top_->localPosition(s1));
0216 GlobalPoint LastStripCenterPointInGlobal = RPCSurface.toGlobal(top_->localPosition(sLast));
0217 double rpcphiFirst = FirstStripCenterPointInGlobal.barePhi();
0218 double rpcphiLast = LastStripCenterPointInGlobal.barePhi();
0219
0220
0221 double diff = rpcphiLast - rpcphiFirst;
0222 double rollphi = (rpcphiFirst + rpcphiLast) * 0.5 * 180. / 3.141592;
0223 double orientation = diff / fabs(diff);
0224 int seg = rpcsrv.segment();
0225
0226
0227
0228
0229 detailstream << std::setw(45) << name << " || phi=" << rollphi << " orientation=" << orientation
0230 << " seg=" << seg << " phi first strip=" << rpcphiFirst << " phi last strip=" << rpcphiLast
0231 << " || " << std::endl;
0232 detailstream << std::setw(45) << name << " || glob (X,Y,Z) first strip =("
0233 << FirstStripCenterPointInGlobal.x() << ", " << FirstStripCenterPointInGlobal.y() << ","
0234 << FirstStripCenterPointInGlobal.z() << ")[cm]"
0235 << " glob (X,Y,Z) last strip =(" << LastStripCenterPointInGlobal.x() << ", "
0236 << LastStripCenterPointInGlobal.y() << "," << LastStripCenterPointInGlobal.z() << ")[cm]"
0237 << std::endl;
0238
0239
0240 ++counterRollsBarrel;
0241 if (rpcId.station() == 4) {
0242 ++counterRollsMB4;
0243 switch (rpcId.ring()) {
0244 case -2:
0245 ++RB4Wm2;
0246 break;
0247 case -1:
0248 ++RB4Wm1;
0249 break;
0250 case 0:
0251 ++RB4W0;
0252 break;
0253 case 1:
0254 ++RB4W1;
0255 break;
0256 case 2:
0257 ++RB4W2;
0258 break;
0259 }
0260 } else
0261 ++counterRollsMB1MB2MB3;
0262 }
0263
0264
0265
0266
0267
0268 else {
0269 const TrapezoidalStripTopology* top_ = dynamic_cast<const TrapezoidalStripTopology*>(&((*r)->topology()));
0270 float s1 = static_cast<float>(1) - 0.5;
0271 float sLast = static_cast<float>(stripsinthisroll) - 0.5;
0272 float stripl = top_->stripLength();
0273 float stripw = top_->pitch();
0274 areaendcap = areaendcap + stripw * stripl * stripsinthisroll;
0275 sumstripwendcap = sumstripwendcap + stripw * stripsinthisroll;
0276 std::string name = rpcsrv.name();
0277
0278
0279
0280 detidsstream << rpcId.rawId() << " \t " << std::setw(24) << name << " \t " << rpcId << std::endl;
0281 detailstream << "All Info " << rpcId.rawId() << " = " << std::setw(24) << name << " || ";
0282 detailstream << " strips: length =" << stripl << "[cm] & width =" << stripw
0283 << "[cm] & number =" << stripsinthisroll;
0284 detailstream << " || area roll =" << stripl * stripw << "[cm^2] || area total endcap = " << areaendcap
0285 << "[cm^2]" << std::endl;
0286
0287
0288 const BoundPlane& RPCSurface = (*r)->surface();
0289 GlobalPoint FirstStripCenterPointInGlobal = RPCSurface.toGlobal(top_->localPosition(s1));
0290 GlobalPoint LastStripCenterPointInGlobal = RPCSurface.toGlobal(top_->localPosition(sLast));
0291 double rpcphiFirst = FirstStripCenterPointInGlobal.barePhi();
0292 double rpcphiLast = LastStripCenterPointInGlobal.barePhi();
0293
0294
0295 double diff = rpcphiLast - rpcphiFirst;
0296 double rollphi = (rpcphiFirst + rpcphiLast) * 0.5 * 180. / 3.141592;
0297 double orientation = diff / fabs(diff);
0298 int seg = rpcsrv.segment();
0299
0300 if (seg == 19)
0301 orientation = orientation * -1;
0302
0303
0304
0305
0306 detailstream << std::setw(45) << name << " || phi=" << rollphi << " orientation=" << orientation
0307 << " seg=" << seg << " phi first strip=" << rpcphiFirst << " phi last strip=" << rpcphiLast
0308 << " || " << std::endl;
0309 detailstream << std::setw(45) << name << " || glob (X,Y,Z) first strip =("
0310 << FirstStripCenterPointInGlobal.x() << ", " << FirstStripCenterPointInGlobal.y() << ","
0311 << FirstStripCenterPointInGlobal.z() << ")[cm]"
0312 << " glob (X,Y,Z) last strip =(" << LastStripCenterPointInGlobal.x() << ", "
0313 << LastStripCenterPointInGlobal.y() << "," << LastStripCenterPointInGlobal.z() << ")[cm]";
0314
0315
0316
0317
0318
0319 bool orientation_ok = false;
0320 if (rpcId.station() == 1) {
0321 if (rpcId.ring() == 1 && orientation * rpcId.region() == 1.0) {
0322 orientation_ok = true;
0323 }
0324 if (rpcId.ring() == 2 && seg % 2 != 0 && orientation * rpcId.region() == 1.0) {
0325 orientation_ok = true;
0326 }
0327 if (rpcId.ring() == 2 && seg % 2 == 0 && orientation * rpcId.region() == -1.0) {
0328 orientation_ok = true;
0329 }
0330 if (rpcId.ring() == 3 && orientation * rpcId.region() == 1.0) {
0331 orientation_ok = true;
0332 }
0333 }
0334 if (rpcId.station() == 2 || rpcId.station() == 4) {
0335 if (orientation * rpcId.region() == -1.0) {
0336 orientation_ok = true;
0337 }
0338 }
0339 if (rpcId.station() == 3) {
0340 if (orientation * rpcId.region() == 1.0) {
0341 orientation_ok = true;
0342 }
0343 }
0344 if (orientation_ok)
0345 detailstream << " OK" << std::endl;
0346 else
0347 detailstream << " WRONG!!!" << std::endl;
0348
0349
0350 ++counterRollsEndCap;
0351 if (rpcId.region() == 1) {
0352 switch (rpcId.station()) {
0353 case 1:
0354 switch (rpcId.ring()) {
0355 case 1:
0356 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0357 break;
0358 case 2:
0359 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0360 break;
0361 case 3:
0362 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0363 break;
0364 }
0365 break;
0366 case 2:
0367 switch (rpcId.ring()) {
0368 case 1:
0369 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0370 break;
0371 case 2:
0372 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0373 break;
0374 case 3:
0375 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0376 break;
0377 }
0378 break;
0379 case 3:
0380 switch (rpcId.ring()) {
0381 case 1:
0382 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0383 break;
0384 case 2:
0385 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0386 break;
0387 case 3:
0388 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0389 break;
0390 }
0391 break;
0392 case 4:
0393 switch (rpcId.ring()) {
0394 case 1:
0395 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0396 break;
0397 case 2:
0398 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0399 break;
0400 case 3:
0401 ++ENDCAProll[rpcId.station()][rpcId.ring()];
0402 break;
0403 }
0404 break;
0405 }
0406 }
0407 }
0408
0409
0410
0411
0412 if (rpcId.region() == 1 && rpcId.station() == 2 &&
0413 (rpcId.sector() == 3 || rpcId.sector() == 4 || rpcId.sector() == 5)) {
0414 ++rollsNearDiskp2;
0415 }
0416 if (rpcId.region() == 1 && rpcId.station() == 3 &&
0417 (rpcId.sector() == 3 || rpcId.sector() == 4 || rpcId.sector() == 2)) {
0418 ++rollsNearDiskp3;
0419 }
0420 for (int strip = 1; strip <= stripsinthisroll; ++strip) {
0421
0422
0423
0424
0425 ++StripsInCMS;
0426 if (rpcId.region() == 0)
0427 ++counterstripsBarrel;
0428 else
0429 ++counterstripsEndCap;
0430 }
0431 }
0432 }
0433 }
0434
0435 detailstream << "Total Number of Strips in CMS=" << StripsInCMS << std::endl;
0436 detailstream << "Total Number of Rolls in CMS=" << RollsInCMS << std::endl;
0437
0438 detailstream << "\n Total Number of Rolls in EndCap= " << counterRollsEndCap << std::endl;
0439 detailstream << "Total Number of Rolls in Barrel= " << counterRollsBarrel << std::endl;
0440
0441 detailstream << "\n Total Number of Strips in Barrel= " << counterstripsBarrel << std::endl;
0442 detailstream << "Total Number of Strips in EndCap= " << counterstripsEndCap << std::endl;
0443
0444 detailstream << "\n Total Number of Rolls in MB4= " << counterRollsMB4 << std::endl;
0445 detailstream << "Total Number of Rolls in MB1,MB2,MB3= " << counterRollsMB1MB2MB3 << std::endl;
0446
0447 detailstream << "RB4 in the barrel:"
0448 << "\n Wheel -2 = " << RB4Wm2 << "\n Wheel -1 = " << RB4Wm1 << "\n Wheel 0 = " << RB4W0
0449 << "\n Wheel +1 = " << RB4W1 << "\n Wheel +2 = " << RB4W2 << std::endl;
0450
0451 detailstream << "In the Endcaps we have the follow distribution of chambers:" << std::endl;
0452 for (int i = 1; i < 5; i++) {
0453 for (int j = 1; j < 4; j++) {
0454 detailstream << " Station " << i << " Ring " << j << " " << ENDCAP[i][j] << " Chambers" << std::endl;
0455 }
0456 }
0457
0458 detailstream << "In the Endcaps we have the follow distribution of rolls:" << std::endl;
0459 for (int i = 1; i < 5; i++) {
0460 for (int j = 1; j < 4; j++) {
0461 detailstream << " Station " << i << " Ring " << j << " " << ENDCAProll[i][j] << " Roll" << std::endl;
0462 }
0463 }
0464
0465 detailstream << "Rolls in Near Disk 2= " << rollsNearDiskp2 << std::endl;
0466 detailstream << "Rolls in Near Disk 3= " << rollsNearDiskp3 << std::endl;
0467
0468 detailstream << "Average Strip Width in Barrel= " << sumstripwbarrel / counterstripsBarrel << std::endl;
0469 detailstream << "Average Strip Width in EndCap= " << sumstripwendcap / counterstripsEndCap << std::endl;
0470
0471 detailstream << "Expected RMS Barrel= " << (sumstripwbarrel / counterstripsBarrel) / sqrt(12) << std::endl;
0472 detailstream << "Expected RMS EndCap= " << (sumstripwendcap / counterstripsEndCap) / sqrt(12) << std::endl;
0473
0474 detailstream << "Area Covered in Barrel = " << areabarrel / 10000 << "m^2" << std::endl;
0475 detailstream << "Area Covered in EndCap = " << areaendcap / 10000 << "m^2" << std::endl;
0476 }
0477
0478
0479 DEFINE_FWK_MODULE(RPCGEO);