File indexing completed on 2024-12-20 03:14:06
0001 #include "Geometry/RPCGeometry/interface/RPCGeometry.h"
0002 #include "Geometry/DTGeometry/interface/DTGeometry.h"
0003 #include "DataFormats/DTRecHit/interface/DTRecSegment4DCollection.h"
0004 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0005 #include "Geometry/Records/interface/MuonGeometryRecord.h"
0006 #include "Geometry/CommonTopologies/interface/RectangularStripTopology.h"
0007 #include "FWCore/Framework/interface/ESHandle.h"
0008 #include "Geometry/RPCGeometry/interface/RPCGeomServ.h"
0009 #include "DataFormats/RPCRecHit/interface/RPCRecHit.h"
0010 #include "DataFormats/RPCRecHit/interface/RPCRecHitCollection.h"
0011 #include "DTSegtoRPC.h"
0012 #include "DTObjectMap.h"
0013 #include "DTStationIndex.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include <ctime>
0017
0018 int distsector(int sector1, int sector2) {
0019 if (sector1 == 13)
0020 sector1 = 4;
0021 if (sector1 == 14)
0022 sector1 = 10;
0023
0024 if (sector2 == 13)
0025 sector2 = 4;
0026 if (sector2 == 14)
0027 sector2 = 10;
0028
0029 int distance = std::abs(sector1 - sector2);
0030 if (distance > 6)
0031 distance = 12 - distance;
0032 return distance;
0033 }
0034
0035 int distwheel(int wheel1, int wheel2) {
0036 int distance = std::abs(wheel1 - wheel2);
0037 return distance;
0038 }
0039
0040 DTSegtoRPC::DTSegtoRPC(edm::ConsumesCollector iC, const edm::ParameterSet& iConfig)
0041 : rpcGeoToken_(iC.esConsumes()), dtGeoToken_(iC.esConsumes()), dtMapToken_(iC.esConsumes()) {
0042 minPhiBX = iConfig.getParameter<int>("minBX");
0043 maxPhiBX = iConfig.getParameter<int>("maxBX");
0044 incldt = true;
0045 incldtMB4 = true;
0046
0047
0048 MinCosAng = 0.85;
0049 MaxD = 80.;
0050 MaxDrb4 = 150.;
0051 MaxDistanceBetweenSegments = 150;
0052 }
0053
0054 std::unique_ptr<RPCRecHitCollection> DTSegtoRPC::thePoints(const DTRecSegment4DCollection* all4DSegments,
0055 const edm::EventSetup& iSetup,
0056 bool debug,
0057 double eyr) {
0058 auto _ThePoints = std::make_unique<RPCRecHitCollection>();
0059 edm::OwnVector<RPCRecHit> RPCPointVector;
0060 std::vector<uint32_t> extrapolatedRolls;
0061
0062 if (all4DSegments->size() > 8) {
0063 if (debug)
0064 LogDebug("DTSegtoRPC") << "Too many segments in this event we are not doing the extrapolation" << std::endl;
0065 } else {
0066 edm::ESHandle<RPCGeometry> rpcGeo = iSetup.getHandle(rpcGeoToken_);
0067 edm::ESHandle<DTGeometry> dtGeo = iSetup.getHandle(dtGeoToken_);
0068 edm::ESHandle<DTObjectMap> dtMap = iSetup.getHandle(dtMapToken_);
0069
0070 std::map<DTChamberId, int> DTSegmentCounter;
0071 DTRecSegment4DCollection::const_iterator segment;
0072
0073 for (segment = all4DSegments->begin(); segment != all4DSegments->end(); ++segment) {
0074 DTSegmentCounter[segment->chamberId()]++;
0075 }
0076
0077 const float bunchCrossTimeDiff = 25.;
0078
0079 if (incldt) {
0080 for (segment = all4DSegments->begin(); segment != all4DSegments->end(); ++segment) {
0081 DTChamberId DTId = segment->chamberId();
0082
0083 if (debug)
0084 LogDebug("DTSegtoRPC") << "DT \t \t This Segment is in Chamber id: " << DTId << std::endl;
0085 if (debug)
0086 LogDebug("DTSegtoRPC") << "DT \t \t Number of segments in this DT = " << DTSegmentCounter[DTId] << std::endl;
0087 if (debug)
0088 LogDebug("DTSegtoRPC") << "DT \t \t Is the only one in this DT? and is not in the 4th Station?" << std::endl;
0089
0090 if (DTSegmentCounter[DTId] != 1 || DTId.station() == 4) {
0091 if (debug)
0092 LogDebug("DTSegtoRPC") << "DT \t \t More than one segment in this chamber, or we are in Station 4"
0093 << std::endl;
0094 continue;
0095 }
0096
0097 int dtWheel = DTId.wheel();
0098 int dtStation = DTId.station();
0099 int dtSector = DTId.sector();
0100
0101 LocalPoint segmentPosition = segment->localPosition();
0102 LocalVector segmentDirection = segment->localDirection();
0103
0104 const GeomDet* gdet = dtGeo->idToDet(segment->geographicalId());
0105 const BoundPlane& DTSurface = gdet->surface();
0106
0107
0108
0109 if (debug)
0110 LogDebug("DTSegtoRPC") << "DT \t \t Is the segment 4D?" << std::endl;
0111
0112 if (segment->dimension() != 4) {
0113 if (debug)
0114 LogDebug("DTSegtoRPC") << "DT \t \t no" << std::endl;
0115 continue;
0116 }
0117
0118 if (debug)
0119 LogDebug("DTSegtoRPC") << "DT \t \t yes" << std::endl;
0120 if (debug)
0121 LogDebug("DTSegtoRPC") << "DT \t \t DT Segment Dimension " << segment->dimension() << std::endl;
0122
0123 float Xo = segmentPosition.x();
0124 float Yo = segmentPosition.y();
0125 float Zo = segmentPosition.z();
0126 float dx = segmentDirection.x();
0127 float dy = segmentDirection.y();
0128 float dz = segmentDirection.z();
0129
0130 float myPhiTime = -9999.;
0131 float myPhiTimeErr = -9999.;
0132 int myPhiBx = -99;
0133
0134 if (segment->hasPhi()) {
0135 myPhiTime = segment->phiSegment()->t0();
0136 myPhiBx = round(segment->phiSegment()->t0() / bunchCrossTimeDiff);
0137 if (debug)
0138 LogDebug("DTSegtoRPC") << "segment t0 = " << myPhiTime << "\tround = " << myPhiBx << std::endl;
0139 }
0140 if (!(myPhiBx <= maxPhiBX && myPhiBx >= minPhiBX))
0141 continue;
0142
0143 if (debug)
0144 LogDebug("DTSegtoRPC") << "Creating the DTIndex" << std::endl;
0145 DTStationIndex theindex(0, dtWheel, dtSector, dtStation);
0146 if (debug)
0147 LogDebug("DTSegtoRPC") << "Getting the Rolls for the given index" << std::endl;
0148 std::set<RPCDetId> rollsForThisDT = dtMap->getRolls(theindex);
0149
0150 if (debug)
0151 LogDebug("DTSegtoRPC") << "DT \t \t Number of rolls for this DT = " << rollsForThisDT.size() << std::endl;
0152
0153 assert(!rollsForThisDT.empty());
0154
0155 if (debug)
0156 LogDebug("DTSegtoRPC") << "DT \t \t Loop over all the rolls asociated to this DT" << std::endl;
0157 for (std::set<RPCDetId>::iterator iteraRoll = rollsForThisDT.begin(); iteraRoll != rollsForThisDT.end();
0158 iteraRoll++) {
0159 const RPCRoll* rollasociated = rpcGeo->roll(*iteraRoll);
0160 RPCDetId rpcId = rollasociated->id();
0161 const BoundPlane& RPCSurface = rollasociated->surface();
0162
0163 RPCGeomServ rpcsrv(rpcId);
0164 std::string nameRoll = rpcsrv.name();
0165
0166 if (debug)
0167 LogDebug("DTSegtoRPC") << "DT \t \t \t RollName: " << nameRoll << std::endl;
0168 if (debug)
0169 LogDebug("DTSegtoRPC") << "DT \t \t \t Doing the extrapolation to this roll" << std::endl;
0170 if (debug)
0171 LogDebug("DTSegtoRPC") << "DT \t \t \t DT Segment Direction in DTLocal " << segmentDirection << std::endl;
0172 if (debug)
0173 LogDebug("DTSegtoRPC") << "DT \t \t \t DT Segment Point in DTLocal " << segmentPosition << std::endl;
0174
0175 GlobalPoint CenterPointRollGlobal = RPCSurface.toGlobal(LocalPoint(0, 0, 0));
0176
0177 LocalPoint CenterRollinDTFrame = DTSurface.toLocal(CenterPointRollGlobal);
0178
0179 if (debug)
0180 LogDebug("DTSegtoRPC") << "DT \t \t \t Center (0,0,0) Roll In DTLocal" << CenterRollinDTFrame << std::endl;
0181 if (debug)
0182 LogDebug("DTSegtoRPC") << "DT \t \t \t Center (0,0,0) of the Roll in Global" << CenterPointRollGlobal
0183 << std::endl;
0184
0185 float D = CenterRollinDTFrame.z();
0186
0187 float X = Xo + dx * D / dz;
0188 float Y = Yo + dy * D / dz;
0189 float Z = D;
0190
0191 const RectangularStripTopology* top_ =
0192 dynamic_cast<const RectangularStripTopology*>(&(rollasociated->topology()));
0193 LocalPoint xmin = top_->localPosition(0.);
0194 if (debug)
0195 LogDebug("DTSegtoRPC") << "DT \t \t \t xmin of this Roll " << xmin << "cm" << std::endl;
0196 LocalPoint xmax = top_->localPosition((float)rollasociated->nstrips());
0197 if (debug)
0198 LogDebug("DTSegtoRPC") << "DT \t \t \t xmax of this Roll " << xmax << "cm" << std::endl;
0199 float rsize = fabs(xmax.x() - xmin.x());
0200 if (debug)
0201 LogDebug("DTSegtoRPC") << "DT \t \t \t Roll Size " << rsize << "cm" << std::endl;
0202 float stripl = top_->stripLength();
0203
0204 float stripw = top_->pitch();
0205
0206 if (debug)
0207 LogDebug("DTSegtoRPC") << "DT \t \t \t Strip Lenght " << stripl << "cm" << std::endl;
0208 if (debug)
0209 LogDebug("DTSegtoRPC") << "DT \t \t \t Strip Width " << stripw << "cm" << std::endl;
0210 if (debug)
0211 LogDebug("DTSegtoRPC") << "DT \t \t \t X Predicted in DTLocal= " << X << "cm" << std::endl;
0212 if (debug)
0213 LogDebug("DTSegtoRPC") << "DT \t \t \t Y Predicted in DTLocal= " << Y << "cm" << std::endl;
0214 if (debug)
0215 LogDebug("DTSegtoRPC") << "DT \t \t \t Z Predicted in DTLocal= " << Z << "cm" << std::endl;
0216
0217 float extrapolatedDistance = sqrt((X - Xo) * (X - Xo) + (Y - Yo) * (Y - Yo) + (Z - Zo) * (Z - Zo));
0218
0219 if (debug)
0220 LogDebug("DTSegtoRPC") << "DT \t \t \t Is the distance of extrapolation less than MaxD? ="
0221 << extrapolatedDistance << "cm"
0222 << "MaxD=" << MaxD << "cm" << std::endl;
0223
0224 if (extrapolatedDistance <= MaxD) {
0225 if (debug)
0226 LogDebug("DTSegtoRPC") << "DT \t \t \t yes" << std::endl;
0227 GlobalPoint GlobalPointExtrapolated = DTSurface.toGlobal(LocalPoint(X, Y, Z));
0228 if (debug)
0229 LogDebug("DTSegtoRPC") << "DT \t \t \t Point ExtraPolated in Global" << GlobalPointExtrapolated
0230 << std::endl;
0231 LocalPoint PointExtrapolatedRPCFrame = RPCSurface.toLocal(GlobalPointExtrapolated);
0232
0233 if (debug)
0234 LogDebug("DTSegtoRPC") << "DT \t \t \t Point Extrapolated in RPCLocal" << PointExtrapolatedRPCFrame
0235 << std::endl;
0236 if (debug)
0237 LogDebug("DTSegtoRPC") << "DT \t \t \t Corner of the Roll = (" << rsize * eyr << "," << stripl * eyr
0238 << ")" << std::endl;
0239 if (debug)
0240 LogDebug("DTSegtoRPC") << "DT \t \t \t Info About the Point Extrapolated in X Abs ("
0241 << fabs(PointExtrapolatedRPCFrame.x()) << ","
0242 << fabs(PointExtrapolatedRPCFrame.y()) << ","
0243 << fabs(PointExtrapolatedRPCFrame.z()) << ")" << std::endl;
0244 if (debug)
0245 LogDebug("DTSegtoRPC") << "DT \t \t \t Does the extrapolation go inside this roll?" << std::endl;
0246
0247 if (fabs(PointExtrapolatedRPCFrame.z()) < 1. && fabs(PointExtrapolatedRPCFrame.x()) < rsize * eyr &&
0248 fabs(PointExtrapolatedRPCFrame.y()) < stripl * eyr) {
0249 if (debug)
0250 LogDebug("DTSegtoRPC") << "DT \t \t \t \t yes" << std::endl;
0251 if (debug)
0252 LogDebug("DTSegtoRPC") << "DT \t \t \t \t Creating the RecHit" << std::endl;
0253
0254 RPCRecHit RPCPoint(rpcId, myPhiBx, PointExtrapolatedRPCFrame);
0255 RPCPoint.setTimeAndError(myPhiTime, myPhiTimeErr);
0256
0257 if (debug)
0258 LogDebug("DTSegtoRPC") << "DT \t \t \t \t Clearing the vector" << std::endl;
0259 RPCPointVector.clear();
0260 if (debug)
0261 LogDebug("DTSegtoRPC") << "DT \t \t \t \t Pushing back" << std::endl;
0262 RPCPointVector.push_back(RPCPoint);
0263 if (debug)
0264 LogDebug("DTSegtoRPC") << "DT \t \t \t \t Putting the vector" << std::endl;
0265 _ThePoints->put(rpcId, RPCPointVector.begin(), RPCPointVector.end());
0266
0267 if (debug)
0268 LogDebug("DTSegtoRPC") << "DT \t \t \t \t Filling container with " << nameRoll
0269 << " Point.x=" << PointExtrapolatedRPCFrame.x()
0270 << " Point.y=" << PointExtrapolatedRPCFrame.y()
0271 << " size=" << RPCPointVector.size() << std::endl;
0272
0273 } else {
0274 if (debug)
0275 LogDebug("DTSegtoRPC") << "DT \t \t \t \t No the prediction is outside of this roll" << std::endl;
0276 }
0277 } else {
0278 if (debug)
0279 LogDebug("DTSegtoRPC") << "DT \t \t \t No, Exrtrapolation too long!, canceled" << std::endl;
0280 }
0281 }
0282 }
0283 }
0284
0285 if (incldtMB4) {
0286 if (all4DSegments->size() > 0) {
0287 if (debug)
0288 LogDebug("DTSegtoRPC") << "MB4 \t \t Loop Over all4DSegments " << all4DSegments->size() << std::endl;
0289 extrapolatedRolls.clear();
0290 for (segment = all4DSegments->begin(); segment != all4DSegments->end(); ++segment) {
0291 DTChamberId DTId = segment->chamberId();
0292
0293 float myPhiTime = -9999.;
0294 float myPhiTimeErr = -9999.;
0295 int myPhiBx = -99;
0296
0297 if (segment->hasPhi()) {
0298 myPhiTime = segment->phiSegment()->t0();
0299 myPhiBx = round(segment->phiSegment()->t0() / bunchCrossTimeDiff);
0300 if (debug)
0301 LogDebug("DTSegtoRPC") << "segment t0 = " << myPhiTime << "\tround = " << myPhiBx << std::endl;
0302 }
0303 if (!(myPhiBx <= maxPhiBX && myPhiBx >= minPhiBX))
0304 continue;
0305
0306 if (debug)
0307 LogDebug("DTSegtoRPC") << "MB4 \t \t This Segment is in Chamber id: " << DTId << std::endl;
0308 if (debug)
0309 LogDebug("DTSegtoRPC") << "MB4 \t \t Number of segments in this DT = " << DTSegmentCounter[DTId]
0310 << std::endl;
0311 if (debug)
0312 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Is the only one in this DT? and is in the Station 4?" << std::endl;
0313
0314 if (DTSegmentCounter[DTId] == 1 && DTId.station() == 4) {
0315 if (debug)
0316 LogDebug("DTSegtoRPC") << "MB4 \t \t \t yes" << std::endl;
0317 int dtWheel = DTId.wheel();
0318 int dtStation = DTId.station();
0319 int dtSector = DTId.sector();
0320
0321 LocalPoint segmentPosition = segment->localPosition();
0322 LocalVector segmentDirection = segment->localDirection();
0323
0324 if (debug)
0325 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t The Segment in MB4 is 2D?" << std::endl;
0326 if (segment->dimension() == 2) {
0327 if (debug)
0328 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t yes" << std::endl;
0329 const LocalVector& segmentDirectionMB4 = segmentDirection;
0330 const LocalPoint& segmentPositionMB4 = segmentPosition;
0331
0332 const BoundPlane& DTSurface4 = dtGeo->idToDet(DTId)->surface();
0333
0334 DTRecSegment4DCollection::const_iterator segMB3;
0335
0336 if (debug)
0337 LogDebug("DTSegtoRPC")
0338 << "MB4 \t \t \t \t Loop on segments in =sector && MB3 && adjacent sectors && y dim=4" << std::endl;
0339 for (segMB3 = all4DSegments->begin(); segMB3 != all4DSegments->end(); ++segMB3) {
0340 DTChamberId dtid3 = segMB3->chamberId();
0341
0342 if (debug)
0343 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Segment in Chamber =" << dtid3 << std::endl;
0344
0345 if (distsector(dtid3.sector(), DTId.sector()) <=
0346 1
0347 && distwheel(dtid3.wheel(), DTId.wheel()) <=
0348 1
0349 && dtid3.station() == 3 && DTSegmentCounter[dtid3] == 1 && segMB3->dimension() == 4) {
0350 if (debug)
0351 LogDebug("DTSegtoRPC")
0352 << "MB4 \t \t \t \t distsector =" << distsector(dtid3.sector(), DTId.sector()) << std::endl;
0353 if (debug)
0354 LogDebug("DTSegtoRPC")
0355 << "MB4 \t \t \t \t distwheel =" << distwheel(dtid3.wheel(), DTId.wheel()) << std::endl;
0356
0357 const GeomDet* gdet3 = dtGeo->idToDet(segMB3->geographicalId());
0358 const BoundPlane& DTSurface3 = gdet3->surface();
0359
0360 LocalVector segmentDirectionMB3 = segMB3->localDirection();
0361 GlobalPoint segmentPositionMB3inGlobal = DTSurface3.toGlobal(segMB3->localPosition());
0362 GlobalPoint segmentPositionMB4inGlobal = DTSurface4.toGlobal(segmentPosition);
0363
0364 LocalVector segDirMB3inMB4Frame = DTSurface4.toLocal(DTSurface3.toGlobal(segmentDirectionMB3));
0365
0366 GlobalVector segDirMB4inGlobalFrame = DTSurface4.toGlobal(segmentDirectionMB4);
0367 GlobalVector segDirMB3inGlobalFrame = DTSurface3.toGlobal(segmentDirectionMB3);
0368
0369 float dx = segDirMB4inGlobalFrame.x();
0370 float dy = segDirMB4inGlobalFrame.y();
0371
0372 float dx3 = segDirMB3inGlobalFrame.x();
0373 float dy3 = segDirMB3inGlobalFrame.y();
0374
0375 double cosAng = fabs(dx * dx3 + dy * dy3 / sqrt((dx3 * dx3 + dy3 * dy3) * (dx * dx + dy * dy)));
0376
0377 if (debug)
0378 LogDebug("DTSegtoRPC")
0379 << "MB4 \t \t \t \t cosAng" << cosAng << "Beetween " << dtid3 << " and " << DTId << std::endl;
0380
0381 if (debug) {
0382 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t dx=" << dx << " dy=" << dy << std::endl;
0383 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t dx3=" << dx3 << " dy3=" << dy << std::endl;
0384 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t cosAng=" << cosAng << std::endl;
0385 }
0386
0387 float DistanceBetweenSegments = ((segmentPositionMB3inGlobal) - (segmentPositionMB4inGlobal)).mag();
0388
0389 if (cosAng > MinCosAng && DistanceBetweenSegments < MaxDistanceBetweenSegments) {
0390 if (debug)
0391 LogDebug("DTSegtoRPC")
0392 << "MB4 \t \t \t \t Distance between segments=" << DistanceBetweenSegments << std::endl;
0393
0394 if (debug)
0395 LogDebug("DTSegtoRPC")
0396 << "MB4 \t \t We found compatible Segments (similar direction and close enough) in " << dtid3
0397 << " and " << DTId << std::endl;
0398
0399 if (dtSector == 13) {
0400 dtSector = 4;
0401 }
0402 if (dtSector == 14) {
0403 dtSector = 10;
0404 }
0405
0406 if (debug)
0407 LogDebug("DTSegtoRPC") << "Creating the DTIndex" << std::endl;
0408 DTStationIndex theindex(0, dtWheel, dtSector, dtStation);
0409 if (debug)
0410 LogDebug("DTSegtoRPC") << "Getting the Rolls for the given index" << std::endl;
0411 std::set<RPCDetId> rollsForThisDT = dtMap->getRolls(theindex);
0412
0413 if (debug)
0414 LogDebug("DTSegtoRPC")
0415 << "MB4 \t \t Number of rolls for this DT = " << rollsForThisDT.size() << std::endl;
0416
0417 assert(!rollsForThisDT.empty());
0418
0419 if (debug)
0420 LogDebug("DTSegtoRPC") << "MB4 \t \t Loop over all the rolls asociated to this DT" << std::endl;
0421 for (std::set<RPCDetId>::iterator iteraRoll = rollsForThisDT.begin();
0422 iteraRoll != rollsForThisDT.end();
0423 iteraRoll++) {
0424 const RPCRoll* rollasociated = rpcGeo->roll(*iteraRoll);
0425 RPCDetId rpcId = rollasociated->id();
0426 const BoundPlane& RPCSurfaceRB4 = rollasociated->surface();
0427
0428 RPCGeomServ rpcsrv(rpcId);
0429 std::string nameRoll = rpcsrv.name();
0430
0431 if (debug)
0432 LogDebug("DTSegtoRPC") << "MB4 \t \t \t RollName: " << nameRoll << std::endl;
0433 if (debug)
0434 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Doing the extrapolation to this roll" << std::endl;
0435
0436 GlobalPoint CenterPointRollGlobal = RPCSurfaceRB4.toGlobal(LocalPoint(0, 0, 0));
0437 LocalPoint CenterRollinMB4Frame = DTSurface4.toLocal(CenterPointRollGlobal);
0438 LocalPoint segmentPositionMB3inMB4Frame =
0439 DTSurface4.toLocal(segmentPositionMB3inGlobal);
0440
0441 LocalVector segmentDirectionMB3inMB4Frame = DTSurface4.toLocal(segDirMB3inGlobalFrame);
0442
0443
0444 float Dxz = CenterRollinMB4Frame.z();
0445 float Xo4 = segmentPositionMB4.x();
0446 float dxl = segmentDirectionMB4.x();
0447 float dzl = segmentDirectionMB4.z();
0448
0449 float X = Xo4 + dxl * Dxz / dzl;
0450 float Z = Dxz;
0451
0452
0453 float Yo34 = segmentPositionMB3inMB4Frame.y();
0454 float dy34 = segmentDirectionMB3inMB4Frame.y();
0455 float dz34 = segmentDirectionMB3inMB4Frame.z();
0456 float Dy =
0457 Dxz -
0458 (segmentPositionMB3inMB4Frame.z());
0459
0460 if (debug)
0461 LogDebug("DTSegtoRPC")
0462 << "MB4 \t \t \t The distance to extrapolate in Y from MB3 is " << Dy << "cm" << std::endl;
0463
0464 float Y = Yo34 + dy34 * Dy / dz34;
0465
0466 const RectangularStripTopology* top_ = dynamic_cast<const RectangularStripTopology*>(
0467 &(rollasociated->topology()));
0468 LocalPoint xmin = top_->localPosition(0.);
0469 LocalPoint xmax = top_->localPosition((float)rollasociated->nstrips());
0470 float rsize = fabs(xmax.x() - xmin.x());
0471 float stripl = top_->stripLength();
0472 float stripw = top_->pitch();
0473
0474 if (debug)
0475 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Strip Lenght " << stripl << "cm" << std::endl;
0476 if (debug)
0477 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Strip Width " << stripw << "cm" << std::endl;
0478
0479 if (debug)
0480 LogDebug("DTSegtoRPC") << "MB4 \t \t \t X Predicted in MB4DTLocal= " << X << "cm" << std::endl;
0481 if (debug)
0482 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Y Predicted in MB4DTLocal= " << Y << "cm" << std::endl;
0483 if (debug)
0484 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Z Predicted in MB4DTLocal= " << Z << "cm" << std::endl;
0485
0486 float extrapolatedDistance = sqrt((Y - Yo34) * (Y - Yo34) + Dy * Dy);
0487
0488 if (debug)
0489 LogDebug("DTSegtoRPC")
0490 << "MB4 \t \t \t segmentPositionMB3inMB4Frame" << segmentPositionMB3inMB4Frame << std::endl;
0491 if (debug)
0492 LogDebug("DTSegtoRPC")
0493 << "MB4 \t \t \t segmentPositionMB4inMB4Frame" << segmentPosition << std::endl;
0494
0495 if (debug)
0496 LogDebug("DTSegtoRPC")
0497 << "MB4 \t \t \t segmentDirMB3inMB4Frame" << segDirMB3inMB4Frame << std::endl;
0498 if (debug)
0499 LogDebug("DTSegtoRPC")
0500 << "MB4 \t \t \t segmentDirMB4inMB4Frame" << segmentDirectionMB4 << std::endl;
0501
0502 if (debug)
0503 LogDebug("DTSegtoRPC")
0504 << "MB4 \t \t \t CenterRB4PositioninMB4Frame" << CenterRollinMB4Frame << std::endl;
0505
0506 if (debug)
0507 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Is the extrapolation distance =" << extrapolatedDistance
0508 << "less than " << MaxDrb4 << std::endl;
0509
0510 if (extrapolatedDistance <= MaxDrb4) {
0511 if (debug)
0512 LogDebug("DTSegtoRPC") << "MB4 \t \t \t yes" << std::endl;
0513
0514 GlobalPoint GlobalPointExtrapolated = DTSurface4.toGlobal(LocalPoint(X, Y, Z));
0515
0516 if (debug)
0517 LogDebug("DTSegtoRPC")
0518 << "MB4 \t \t \t Point ExtraPolated in Global" << GlobalPointExtrapolated << std::endl;
0519
0520 LocalPoint PointExtrapolatedRPCFrame = RPCSurfaceRB4.toLocal(GlobalPointExtrapolated);
0521
0522 if (debug)
0523 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Point Extrapolated in RPCLocal"
0524 << PointExtrapolatedRPCFrame << std::endl;
0525 if (debug)
0526 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Corner of the Roll = (" << rsize * eyr << ","
0527 << stripl * eyr << ")" << std::endl;
0528 if (debug)
0529 LogDebug("DTSegtoRPC")
0530 << "MB4 \t \t \t Info About the Point Extrapolated in X Abs ("
0531 << fabs(PointExtrapolatedRPCFrame.x()) << "," << fabs(PointExtrapolatedRPCFrame.y())
0532 << "," << fabs(PointExtrapolatedRPCFrame.z()) << ")" << std::endl;
0533
0534 if (debug)
0535 LogDebug("DTSegtoRPC")
0536 << "MB4 \t \t \t Does the extrapolation go inside this roll?" << std::endl;
0537
0538 if (fabs(PointExtrapolatedRPCFrame.z()) < 5. &&
0539 fabs(PointExtrapolatedRPCFrame.x()) < rsize * eyr &&
0540 fabs(PointExtrapolatedRPCFrame.y()) < stripl * eyr) {
0541 if (debug)
0542 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t yes" << std::endl;
0543 if (debug)
0544 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Creating the RecHit" << std::endl;
0545
0546 RPCRecHit RPCPointMB4(rpcId, myPhiBx, PointExtrapolatedRPCFrame);
0547 RPCPointMB4.setTimeAndError(myPhiTime, myPhiTimeErr);
0548
0549 if (debug)
0550 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Clearing the RPCPointVector" << std::endl;
0551 RPCPointVector.clear();
0552 if (debug)
0553 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Pushing Back" << std::endl;
0554 RPCPointVector.push_back(RPCPointMB4);
0555 if (debug)
0556 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Putting for " << rpcId << std::endl;
0557 if (debug)
0558 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Filling container with " << nameRoll
0559 << " Point.x=" << PointExtrapolatedRPCFrame.x()
0560 << " Point.y=" << PointExtrapolatedRPCFrame.y()
0561 << " size=" << RPCPointVector.size() << std::endl;
0562 if (debug)
0563 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t Number of rolls already extrapolated in RB4 = "
0564 << extrapolatedRolls.size() << std::endl;
0565 if (find(extrapolatedRolls.begin(), extrapolatedRolls.end(), rpcId.rawId()) ==
0566 extrapolatedRolls.end()) {
0567 extrapolatedRolls.push_back(rpcId.rawId());
0568 _ThePoints->put(rpcId, RPCPointVector.begin(), RPCPointVector.end());
0569 } else {
0570 if (debug)
0571 LogDebug("DTSegtoRPC")
0572 << "MB4 \t \t \t \t roll already extrapolated " << rpcId << std::endl;
0573 }
0574 if (debug)
0575 LogDebug("DTSegtoRPC")
0576 << "MB4 \t \t \t \t Extrapolations done after this point = " << extrapolatedRolls.size()
0577 << std::endl;
0578 if (debug)
0579 for (uint32_t m = 0; m < extrapolatedRolls.size(); m++)
0580 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t" << extrapolatedRolls.at(m) << std::endl;
0581 } else {
0582 if (debug)
0583 LogDebug("DTSegtoRPC")
0584 << "MB4 \t \t \t \t No the prediction is outside of this roll" << std::endl;
0585 }
0586 }
0587 else {
0588 if (debug)
0589 LogDebug("DTSegtoRPC") << "MB4 \t \t \t No, Exrtrapolation too long!, canceled" << std::endl;
0590 }
0591 }
0592 } else {
0593 if (debug)
0594 LogDebug("DTSegtoRPC")
0595 << "MB4 \t \t \t \t I found segments in MB4 and MB3 adjacent wheel and/or sector but "
0596 "not compatibles, Diferent Directions"
0597 << std::endl;
0598 }
0599 } else {
0600 if (debug)
0601 LogDebug("DTSegtoRPC")
0602 << "MB4 \t \t \t No the same station or same wheel or segment dim in mb3 not 4D" << std::endl;
0603 }
0604 }
0605 } else {
0606 if (debug)
0607 LogDebug("DTSegtoRPC") << "MB4 \t \t \t Is NOT a 2D Segment" << std::endl;
0608 }
0609 } else {
0610 if (debug)
0611 LogDebug("DTSegtoRPC") << "MB4 \t \t \t \t There is not just one segment or is not in station 4"
0612 << std::endl;
0613 }
0614 }
0615 } else {
0616 if (debug)
0617 LogDebug("DTSegtoRPC") << "MB4 \t This event doesn't have 4D Segment" << std::endl;
0618 }
0619 }
0620 }
0621
0622 return _ThePoints;
0623 }