#include <DQM/RPCMonitorDigi/interface/RPCMonitorEfficiency.h>
Public Member Functions | |
virtual void | analyze (const edm::Event &, const edm::EventSetup &) |
std::map< std::string, MonitorElement * > | bookDetUnitMEEff (RPCDetId &detId) |
Booking of MonitoringElemnt for one RPCDetId (= roll). | |
virtual void | endJob (void) |
RPCMonitorEfficiency (const edm::ParameterSet &) | |
Constructor. | |
~RPCMonitorEfficiency () | |
Destructor. | |
Private Attributes | |
std::vector< uint32_t > | _idList |
std::vector< std::map < RPCDetId, int > > | counter |
DQMStore * | dbe |
back-end interface | |
bool | debug |
std::string | digiLabel |
std::string | EffRootFileName |
bool | EffSaveRootFile |
int | EffSaveRootFileEventsInterval |
MonitorElement * | h1 |
TH1F * | hPositionX |
std::map< uint32_t, std::map < std::string, MonitorElement * > > | meCollection |
std::string | nameInLog |
std::ofstream | ofrej |
std::string | theRecHits4DLabel |
std::vector< int > | totalcounter |
Definition at line 34 of file RPCMonitorEfficiency.h.
RPCMonitorEfficiency::RPCMonitorEfficiency | ( | const edm::ParameterSet & | pset | ) | [explicit] |
Constructor.
get hold of back-end interface
Definition at line 71 of file RPCMonitorEfficiency.cc.
References _idList, counter, dbe, digiLabel, EffRootFileName, EffSaveRootFile, EffSaveRootFileEventsInterval, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), theRecHits4DLabel, and totalcounter.
00071 { 00072 std::map<RPCDetId, int> buff; 00073 counter.clear(); 00074 counter.reserve(3); 00075 //std::cout <<" Buff 1"<<std::endl; 00076 counter.push_back(buff); 00077 //std::cout <<" Buff 2"<<std::endl; 00078 counter.push_back(buff); 00079 //std::cout <<" Buff 3"<<std::endl; 00080 counter.push_back(buff); 00081 totalcounter.clear(); 00082 totalcounter.reserve(3); 00083 totalcounter[0]=0; 00084 totalcounter[1]=0; 00085 totalcounter[2]=0; 00086 theRecHits4DLabel = pset.getParameter<std::string>("recHits4DLabel"); 00087 digiLabel=pset.getParameter<std::string>("digiLabel"); 00088 EffSaveRootFile = pset.getUntrackedParameter<bool>("EffSaveRootFile", false); 00089 EffSaveRootFileEventsInterval = pset.getUntrackedParameter<int>("EffEventsInterval", 10000); 00090 EffRootFileName = pset.getUntrackedParameter<std::string>("EffRootFileName", "RPCEfficiency.root"); 00091 00093 dbe = edm::Service<DQMStore>().operator->(); 00094 00095 00096 dbe->showDirStructure(); 00097 00098 _idList.clear(); 00099 //ofrej.open("rejected.txt"); 00100 }
RPCMonitorEfficiency::~RPCMonitorEfficiency | ( | ) |
void RPCMonitorEfficiency::analyze | ( | const edm::Event & | iEvent, | |
const edm::EventSetup & | iSetup | |||
) | [virtual] |
Implements edm::EDAnalyzer.
Definition at line 102 of file RPCMonitorEfficiency.cc.
References _idList, angle(), b3, bookDetUnitMEEff(), BoundSurface::bounds(), counter, GenMuonPlsPt100GeV_cfg::cout, funct::D, digiLabel, lat::endl(), edm::EventSetup::get(), edm::Event::getByLabel(), edm::Event::id(), RPCRoll::id(), it, RPCDetId::layer(), RectangularStripTopology::localPosition(), meCollection, RPCRoll::nstrips(), DetId::rawId(), RPCDetId::region(), res, RPCDetId::ring(), RPCDetId::roll(), RPCDetId::sector(), DTChamberId::sector(), RPCDetId::station(), DTChamberId::station(), RPCRoll::strip(), RectangularStripTopology::stripLength(), RPCDetId::subsector(), GeomDet::surface(), theRecHits4DLabel, Bounds::thickness(), Surface::toGlobal(), GloballyPositioned< T >::toLocal(), RPCRoll::topology(), totalcounter, w3, w4, muonGeometry::wheel, DTChamberId::wheel(), PV3DBase< T, PVType, FrameType >::x(), X, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().
00102 { 00103 std::map<RPCDetId, int> buff; 00104 char layerLabel[128]; 00105 char meIdRPC [128]; 00106 char meIdDT [128]; 00107 float dx=0.,dy=0.,dz=0.,Xo=0.,Yo=0.,X=0.,Y=0.,Z=0.,p1x=0.,p2x=0.,p3x=0.,p4x=0.,p1z=0.,p2z=0.,p3z=0.,p4z=0.,dx3=0.,dy3=0.,dz3=0.,Xo3=0.,Yo3=0.,x3=0.,x4=0.,z3=0.,z4=0.,m3=0.,m4=0.,xc=0.,zc=0.,b3=0.,b4=0.,w3=0.,w4=0.; 00108 00109 float widestrip=5.; 00110 float widestripsRB4=8.; 00111 float circError=3.; 00112 float angle=0.01; 00113 00114 edm::ESHandle<DTGeometry> dtGeo; 00115 iSetup.get<MuonGeometryRecord>().get(dtGeo); 00116 00117 edm::ESHandle<RPCGeometry> rpcGeo; 00118 iSetup.get<MuonGeometryRecord>().get(rpcGeo); 00119 00120 edm::Handle<DTRecSegment4DCollection> all4DSegments; 00121 iEvent.getByLabel(theRecHits4DLabel, all4DSegments); 00122 00123 edm::Handle<RPCDigiCollection> rpcDigis; 00124 iEvent.getByLabel(digiLabel, rpcDigis); 00125 00126 std::map<DTStationIndex,std::set<RPCDetId> > rollstore; 00127 for (TrackingGeometry::DetContainer::const_iterator it=rpcGeo->dets().begin();it<rpcGeo->dets().end();it++){ 00128 RPCRoll* ir = dynamic_cast<RPCRoll*>(*it); 00129 RPCDetId rpcId = ir->id(); 00130 int region=rpcId.region(); 00131 int wheel=rpcId.ring(); 00132 int sector=rpcId.sector(); 00133 int station=rpcId.station(); 00134 00135 DTStationIndex ind(region,wheel,sector,station); 00136 00137 std::set<RPCDetId> myrolls; 00138 if (rollstore.find(ind)!=rollstore.end()){ 00139 myrolls=rollstore[ind]; 00140 } 00141 myrolls.insert(rpcId); 00142 rollstore[ind]=myrolls; 00143 } 00144 00145 00146 if(all4DSegments->size()>0){ 00147 //std::cout<<"Number of Segments in this event = "<<all4DSegments->size()<<std::endl; 00148 00149 00150 std::map<DTChamberId,int> scounter; 00151 DTRecSegment4DCollection::const_iterator segment; 00152 00153 00154 for (segment = all4DSegments->begin();segment!=all4DSegments->end(); ++segment){ 00155 scounter[segment->chamberId()]++; 00156 } 00157 00158 //std::cout<<"Loop over all the 4D Segments"<<std::endl; 00159 //loop over all the 4D Segments 00160 00161 for (segment = all4DSegments->begin(); segment != all4DSegments->end(); ++segment){ 00162 DTChamberId DTId = segment->chamberId(); 00163 //std::cout<<"\t This Segment is in Chamber id: "<<DTId<<std::endl; 00164 //std::cout<<"\t Number of segments in this DT = "<<scounter[DTId]<<std::endl; 00165 //std::cout<<"\t DT Segment Dimension "<<segment->dimension()<<std::endl; 00166 //std::cout<<"\t Is the only in this DT?"<<std::endl; 00167 00168 //there must be only one segment per Chamber 00169 if(scounter[DTId] == 1){ 00170 //std::cout<<"\t \t yes"<<std::endl; 00171 int dtWheel = DTId.wheel(); 00172 int dtStation = DTId.station(); 00173 int dtSector = DTId.sector(); 00174 00175 LocalPoint segmentPosition= segment->localPosition(); 00176 LocalVector segmentDirection=segment->localDirection(); 00177 00178 const GeomDet* gdet=dtGeo->idToDet(segment->geographicalId()); 00179 const BoundPlane & DTSurface = gdet->surface(); 00180 00181 //check if the dimension of the segment is 4 00182 00183 if(segment->dimension()==4){ 00184 Xo=segmentPosition.x(); 00185 Yo=segmentPosition.y(); 00186 dx=segmentDirection.x(); 00187 dy=segmentDirection.y(); 00188 dz=segmentDirection.z(); 00189 //std::cout<<"\t \t Loop over all the rolls asociated to this DT"<<std::endl; 00190 std::set<RPCDetId> rollsForThisDT = 00191 rollstore[DTStationIndex(0,dtWheel,dtSector,dtStation)]; 00192 //Loop over all the rolls 00193 for (std::set<RPCDetId>::iterator iteraRoll = rollsForThisDT.begin();iteraRoll != rollsForThisDT.end(); iteraRoll++){ 00194 const RPCRoll* rollasociated = rpcGeo->roll(*iteraRoll); 00195 //To get the roll's surface 00196 const BoundPlane & RPCSurface = rollasociated->surface(); 00197 //std::cout<<"\t \t RollID: "<<rollasociated->id()<<std::endl; 00198 //std::cout<<"\t \t Doing the extrapolation"<<std::endl; 00199 //std::cout<<"\t \t DT Segment Direction in DTLocal "<<segmentDirection<<std::endl; 00200 //std::cout<<"\t \t DT Segment Point in DTLocal "<<segmentPosition<<std::endl; 00201 00202 GlobalPoint CenterPointRollGlobal = RPCSurface.toGlobal(LocalPoint(0,0,0)); 00203 //std::cout<<"\t \t Center (0,0,0) of the Roll in Global"<<CenterPointRollGlobal<<std::endl; 00204 00205 LocalPoint CenterRollinDTFrame = DTSurface.toLocal(CenterPointRollGlobal); 00206 //std::cout<<"\t \t Center (0,0,0) Roll In DTLocal"<<CenterRollinDTFrame<<std::endl; 00207 00208 float D=CenterRollinDTFrame.z(); 00209 //std::cout<<"\t \t D="<<D<<"cm"<<std::endl; 00210 00211 X=Xo+dx*D/dz; 00212 Y=Yo+dy*D/dz; 00213 Z=D; 00214 00215 const RectangularStripTopology* top_= dynamic_cast<const RectangularStripTopology*> (&(rollasociated->topology())); 00216 LocalPoint xmin = top_->localPosition(0.); 00217 LocalPoint xmax = top_->localPosition((float)rollasociated->nstrips()); 00218 float rsize = fabs( xmax.x()-xmin.x() )*0.5; 00219 float stripl = top_->stripLength(); 00220 00221 //std::cout<<"\t \t X Predicted in DTLocal= "<<X<<"cm"<<std::endl; 00222 //std::cout<<"\t \t Y Predicted in DTLocal= "<<Y<<"cm"<<std::endl; 00223 //std::cout<<"\t \t Z Predicted in DTLocal= "<<Z<<"cm"<<std::endl; 00224 00225 GlobalPoint GlobalPointExtrapolated = 00226 DTSurface.toGlobal(LocalPoint(X,Y,Z)); 00227 //std::cout<<"\t \t Point ExtraPolated in Global"<<GlobalPointExtrapolated<< std::endl; 00228 00229 LocalPoint PointExtrapolatedRPCFrame = 00230 RPCSurface.toLocal(GlobalPointExtrapolated); 00231 //std::cout<<"\t \t Point Extrapolated in RPCLocal"<<PointExtrapolatedRPCFrame<< std::endl; 00232 //std::cout<<"\t \t Does the extrapolation go inside this roll?"<<std::endl; 00233 //conditions to find the right roll to extrapolate 00234 00235 if(fabs(PointExtrapolatedRPCFrame.z()) < 0.01 && fabs(PointExtrapolatedRPCFrame.x()) < rsize && fabs(PointExtrapolatedRPCFrame.y()) < stripl*0.5){ 00236 //std::cout<<"\t \t \t yes"<<std::endl; 00237 //getting the number of the strip 00238 const float stripPredicted = 00239 rollasociated->strip(LocalPoint(PointExtrapolatedRPCFrame.x(),PointExtrapolatedRPCFrame.y(),0.)); 00240 00241 //std::cout<<"\t \t \t Candidate"<<rollasociated->id()<<" "<<"(from DT Segment) STRIP---> "<<stripPredicted<< std::endl; 00242 00243 //--------- HISTOGRAM STRIP PREDICTED FROM DT ------------------- 00244 00245 RPCDetId rollId = rollasociated->id(); 00246 uint32_t id = rollId.rawId(); 00247 00248 RPCDetId otherRollId1,otherRollId2; 00249 00250 if(rollId.roll() == 1){ 00251 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),2); 00252 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),3); 00253 otherRollId1 = tempRollId1; 00254 otherRollId2 = tempRollId2; 00255 } 00256 else if(rollId.roll() == 2){ 00257 00258 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),1); 00259 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),3); 00260 otherRollId1 = tempRollId1; 00261 otherRollId2 = tempRollId2; 00262 } 00263 else if(rollId.roll() == 3){ 00264 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),1); 00265 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),2); 00266 otherRollId1 = tempRollId1; 00267 otherRollId2 = tempRollId2; 00268 } 00269 00270 RPCDigiCollection::Range rpcRangeDigi1=rpcDigis->get(otherRollId1); 00271 RPCDigiCollection::Range rpcRangeDigi2=rpcDigis->get(otherRollId2); 00272 00273 _idList.push_back(id); 00274 00275 char detUnitLabel[128]; 00276 sprintf(detUnitLabel ,"%d",id); 00277 sprintf(layerLabel ,"layer%d_subsector%d_roll%d",rollId.layer(),rollId.subsector(),rollId.roll()); 00278 00279 std::map<uint32_t, std::map<std::string,MonitorElement*> >::iterator meItr = meCollection.find(id); 00280 if (meItr == meCollection.end()){ 00281 meCollection[id] = bookDetUnitMEEff(rollId); 00282 //std::cout << "\t \t \t Create new histograms for "<<layerLabel<<std::endl; 00283 } 00284 00285 std::map<std::string, MonitorElement*> meMap=meCollection[id]; 00286 sprintf(meIdDT,"ExpectedOccupancyFromDT_%s",detUnitLabel); 00287 meMap[meIdDT]->Fill(stripPredicted); 00288 00289 sprintf(meIdDT,"ExpectedOccupancy2DFromDT_%s",detUnitLabel); 00290 meMap[meIdDT]->Fill(stripPredicted,Y); 00291 00292 //std::cout << "\t \t \t One for counterPREDICT"<<std::endl; 00293 totalcounter[0]++; 00294 buff=counter[0]; 00295 buff[rollId]++; 00296 counter[0]=buff; 00297 //------------------------------------------------------------------- 00298 00299 //std::cout<<"\t \t \t We have a Candidate let's see in the digis!"<<std::endl; 00300 00301 bool anycoincidence=false; 00302 int stripDetected = 0; 00303 RPCDigiCollection::Range rpcRangeDigi=rpcDigis->get(rollasociated->id()); 00304 00305 int stripCounter = 0; 00306 00307 for (RPCDigiCollection::const_iterator digiIt = rpcRangeDigi.first;digiIt!=rpcRangeDigi.second;++digiIt){//loop over the digis in the event 00308 stripCounter++; 00309 //std::cout<<"\t \t \t \t Digi "<<*digiIt<<std::endl;//print the digis in the event 00310 stripDetected=digiIt->strip(); 00311 00312 float res = (float)(stripDetected) - stripPredicted; 00313 sprintf(meIdRPC,"RPCResiduals_%s",detUnitLabel); 00314 meMap[meIdRPC]->Fill(res); 00315 00316 sprintf(meIdRPC,"RPCResiduals2D_%s",detUnitLabel); 00317 meMap[meIdRPC]->Fill(res,Y); 00318 00319 if(res > 7) 00320 std::cout<<"STRANGE "<<"EVENTO NUM = "<<iEvent.id().event()<<" Residuo = "<<res<<" Strip Num = "<<stripDetected<<" Strip totali = "<<stripCounter<<std::endl; 00321 //compare the strip Detected with the predicted 00322 if(fabs((float)(stripDetected) - stripPredicted) < widestrip){ 00323 //std::cout <<"\t \t \t \t COINCEDENCE Predict " 00324 // <<stripPredicted<<" Detect " 00325 // <<stripDetected<<std::endl; 00326 anycoincidence=true; 00327 //break;//funciona solo para hacerlo mas rapido 00328 //We can not divide two diferents things 00329 } 00330 } 00331 if (anycoincidence) { 00332 sprintf(meIdRPC,"ExpectedOccupancyFromDT_forCrT_%s",detUnitLabel); 00333 meMap[meIdRPC]->Fill(stripPredicted); 00334 sprintf(meIdRPC,"RealDetectedOccupancy_%s",detUnitLabel); 00335 meMap[meIdRPC]->Fill(stripDetected); 00336 sprintf(meIdRPC,"YExpectedOccupancyFromDT_%s",detUnitLabel); 00337 meMap[meIdRPC]->Fill(Y); 00338 00339 for (RPCDigiCollection::const_iterator digiIt1 = rpcRangeDigi1.first;digiIt1!=rpcRangeDigi1.second;++digiIt1){ 00340 00341 if(fabs(stripDetected - digiIt1->strip()) <= 1){ 00342 00343 sprintf(meIdRPC,"XCrossTalk_1_%s",detUnitLabel); 00344 meMap[meIdRPC]->Fill(stripPredicted); 00345 00346 sprintf(meIdRPC,"XDetectCrossTalk_1_%s",detUnitLabel); 00347 meMap[meIdRPC]->Fill(stripDetected); 00348 00349 sprintf(meIdRPC,"YCrossTalk_1_%s",detUnitLabel); 00350 meMap[meIdRPC]->Fill(Y); 00351 break; 00352 } 00353 } 00354 00355 for (RPCDigiCollection::const_iterator digiIt2 = rpcRangeDigi2.first;digiIt2!=rpcRangeDigi2.second;++digiIt2){ 00356 00357 if(fabs(stripDetected - digiIt2->strip()) <= 1){ 00358 00359 sprintf(meIdRPC,"XCrossTalk_2_%s",detUnitLabel); 00360 meMap[meIdRPC]->Fill(stripPredicted); 00361 00362 sprintf(meIdRPC,"XDetectCrossTalk_2_%s",detUnitLabel); 00363 meMap[meIdRPC]->Fill(stripDetected); 00364 00365 sprintf(meIdRPC,"YCrossTalk_2_%s",detUnitLabel); 00366 meMap[meIdRPC]->Fill(Y); 00367 break; 00368 } 00369 } 00370 00371 sprintf(meIdRPC,"RPCDataOccupancy_%s",detUnitLabel); 00372 meMap[meIdRPC]->Fill(stripPredicted); 00373 00374 sprintf(meIdRPC,"RPCDataOccupancy2D_%s",detUnitLabel); 00375 meMap[meIdRPC]->Fill(stripPredicted,Y); 00376 00377 totalcounter[1]++; 00378 buff=counter[1]; 00379 buff[rollId]++; 00380 counter[1]=buff; 00381 } 00382 else { 00383 //std::cout <<"\t \t \t \t XXXXX THIS PREDICTION DOESN'T HAVE ANY CORRESPONDENCE WITH THE DATA"<<std::endl; 00384 totalcounter[2]++; 00385 buff=counter[2]; 00386 buff[rollId]++; 00387 counter[2]=buff; 00388 //std::cout << "\t \t \t \t One for counterFAIL"<<std::endl; 00389 //ofrej<<"Wh "<<dtWheel<<"\t St "<<dtStation 00390 // <<"\t Se "<<dtSector<<"\t Event " 00391 // <<iEvent.id().event()<<std::endl; 00392 } 00393 } 00394 else { 00395 //std::cout<<"\t \t \t no"<<std::endl; 00396 }//Condition for the right match 00397 }//loop over all the rolls 00398 00399 00400 // dedicated RB4 analysis part, that misses DT 4D segments 00401 00402 00403 }else if(segment->dimension()==2&&dtStation==4){ 00404 00405 LocalVector segmentDirectionMB4=segmentDirection; 00406 LocalPoint segmentPositionMB4=segmentPosition; 00407 00408 00409 //std::cout<<"\t \t 2D in RB4"<<DTId<<" with D="<<segment->dimension()<<localDirection<<segmentPositionMB4<<std::endl; 00410 bool compatiblesegments=false; 00411 Xo=segmentPositionMB4.x(); 00412 dx=segmentDirectionMB4.x(); 00413 dz=segmentDirectionMB4.z(); 00414 //std::cout<<"\t \t Loop over all the segments"<<std::endl; 00415 DTRecSegment4DCollection::const_iterator segMB3; 00416 00417 const BoundPlane& DTSurface4 = dtGeo->idToDet(DTId)->surface(); 00418 w4 = DTSurface4.bounds().thickness()*0.5; // along local Z 00419 00420 for(segMB3=all4DSegments->begin();segMB3!=all4DSegments->end();++segMB3){ 00421 DTChamberId dtid = segMB3->chamberId(); 00422 if(dtid.station()==3){ 00423 const GeomDet* gdet3=dtGeo->idToDet(segMB3->geographicalId()); 00424 const BoundPlane & DTSurface3 = gdet3->surface(); 00425 w3 = DTSurface3.bounds().thickness()*0.5; // along local Z 00426 00427 dx3=segMB3->localDirection().x(); 00428 dy3=segMB3->localDirection().y(); 00429 dz3=segMB3->localDirection().z(); 00430 00431 //LocalVector segDirMB4inMB3Frame=DTSurface3.toLocal(DTSurface4.toGlobal(segmentDirectionMB4)); 00432 00433 if(fabs(dx-dx3)<=angle&&fabs(dz-dz3)<=angle){//same direction? 00434 compatiblesegments=true; 00435 }else{ 00436 //They don't have the same local dir, and the segments in diferent sectors,compatibles with circle? 00437 Xo3=segMB3->localPosition().x(); 00438 Yo3=segMB3->localPosition().y(); 00439 00440 //Do we have segments compatibles with a muon? 00441 00442 // pos=seg.position() 00443 // dir=seg.direction() 00444 // dest =pos + dir*height*cos(dir.theta()) 00445 00446 //std::cout<<"Ancho de MB3 w3="<<w3<<" |Ancho de MB4 w4= "<<w4<<std::endl; 00447 //std::cout<<"Informacion Inicial MB4"<<segment->localDirection()<<" "<<segment->localPosition()<<std::endl; 00448 //std::cout<<"Informacion Inicial MB3"<<segMB3->localDirection()<<" "<<segMB3->localPosition()<<std::endl; 00449 compatiblesegments=false; 00450 00451 //This in MB3 Local 00452 00453 p1x=Xo3+dx3*w3/dz3; 00454 p1z=w3; 00455 00456 p2x=Xo3-dx3*w3/dz3; 00457 p2z=-w3; 00458 00459 //This in MB4 Local 00460 00461 p3x=Xo+dx*w4/dz; 00462 p3z=w4; 00463 00464 p4x=Xo-dx*w4/dz; 00465 p4z=-w4; 00466 00467 LocalPoint P1=LocalPoint(p1x,0,p1z); 00468 LocalPoint P2=LocalPoint(p2x,0,p2z); 00469 LocalPoint P3=LocalPoint(p3x,0,p3z); 00470 LocalPoint P4=LocalPoint(p4x,0,p4z); 00471 00472 //std::cout<<"Points in MB3="<<P1<<" --- "<<P2<<std::endl; 00473 //std::cout<<"Points in MB4="<<P3<<" --- "<<P4<<std::endl; 00474 00475 //Now we have to convert to global 00476 00477 LocalPoint P1g=P1; 00478 LocalPoint P2g=P2; 00479 00480 LocalPoint P3g=DTSurface3.toLocal(DTSurface.toGlobal(P3)); 00481 LocalPoint P4g=DTSurface3.toLocal(DTSurface.toGlobal(P4)); 00482 00483 //std::cout<<"Points in MB3 ="<<P1g<<" --- "<<P2g<<std::endl; 00484 //std::cout<<"Points in MB4 in MB3Frame="<<P3g<<" --- "<<P4g<<std::endl; 00485 00486 00487 float dx3g=P1g.x()-P2g.x(); 00488 float dz3g=P1g.z()-P2g.z(); 00489 00490 float dxg=P3g.x()-P4g.x(); 00491 float dzg=P3g.z()-P4g.z(); 00492 00493 //std::cout<<"dx3g="<<dx3g<<" dz3g="<<dz3g<<std::endl; 00494 //std::cout<<"dxg="<<dxg<<" dzg="<<dzg<<std::endl; 00495 00496 m3=-dx3g/dz3g; 00497 m4=-dxg/dzg; 00498 00499 x3=(P1g.x()+P2g.x())*0.5; 00500 z3=(P1g.z()+P2g.z())*0.5; 00501 00502 x4=(P3g.x()+P4g.x())*0.5; 00503 z4=(P3g.z()+P4g.z())*0.5; 00504 00505 b3=z3-m3*x3; 00506 b4=z4-m4*x4; 00507 00508 if(m3!=m4){ 00509 xc=(b4-b3)/(m3-m4); 00510 zc=m3*xc+b3; 00511 00512 GlobalPoint Pc=GlobalPoint(xc,0,zc); 00513 00514 //std::cout<<Pc<<std::endl; 00515 00516 float distance=fabs((GlobalPoint(P2g.x(),0,P2g.z())-GlobalPoint(Pc.x(),0,Pc.z())).mag()-(GlobalPoint(P3g.x(),0,P3g.z())-GlobalPoint(Pc.x(),0,Pc.z())).mag()); 00517 00518 if(distance<circError){ 00519 compatiblesegments=true; 00520 //std::cout<<"YES in the same circle... p2 a C="<<P1g<<P2g<<P3g<<P4g<<"Distancia "<<distance<<std::endl; 00521 } 00522 else{ 00523 //std::cout<<"NOT in the same circle... p2 a C="<<P1g<<P2g<<P3g<<P4g<<"Distancia "<<distance<<std::endl; 00524 compatiblesegments=false; 00525 } 00526 } 00527 else{ 00528 //std::cout<<"We have the same slope m3="<<m3<<" m4="<<m4<<std::endl; 00529 if(fabs(b3-b4)<=circError){ 00530 compatiblesegments=true; 00531 //std::cout<<"and the segments are in a line"<<std::endl; 00532 }else{std::cout<<"But we don't have the same intercept b3="<<b3<<" b4="<<b4<<std::endl;} 00533 compatiblesegments=false; 00534 system("sleep 30"); 00535 } 00536 } 00537 00538 //conditions in MB3 00539 if(scounter[dtid]==1 && compatiblesegments){ 00540 //std::cout<<"********\t \t \t In the same event there is a segment in RB3 "<<dtid<<" with D="<<segMB3->dimension()<<segMB3->localDirection()<<segMB3->localPosition()<<"scounter "<<scounter[dtid]<<std::endl; 00541 00542 std::set<RPCDetId> rollsForThisDT = 00543 rollstore[DTStationIndex(0,dtWheel,dtSector,dtStation)]; 00544 //Loop over all the rolls asociated to RB4 00545 for (std::set<RPCDetId>::iterator iteraRoll=rollsForThisDT.begin();iteraRoll != rollsForThisDT.end(); iteraRoll++){ 00546 const RPCRoll* rollasociated = rpcGeo->roll(*iteraRoll); 00547 //To get the roll's surface 00548 const BoundPlane & RPCSurfaceRB4 = rollasociated->surface(); 00549 00550 const GeomDet* gdet=dtGeo->idToDet(segMB3->geographicalId()); 00551 const BoundPlane & DTSurfaceMB3 = gdet->surface(); 00552 00553 //std::cout<<"\t \t \t RollID: should be RB4"<<rollasociated->id()<<std::endl; 00554 //std::cout<<"\t \t \t Making the extrapolation"<<std::endl; 00555 //std::cout<<"\t \t \t DT Segment Direction in MB3 DTLocal "<<segMB3->localDirection()<<std::endl; 00556 //std::cout<<"\t \t \t DT Segment Point in MB3 DTLocal "<<segMB3->localPosition()<<std::endl; 00557 00558 GlobalPoint CenterPointRollGlobal=RPCSurfaceRB4.toGlobal(LocalPoint(0,0,0)); 00559 //std::cout<<"\t \t \t Center (0,0,0) of the RB4 Roll in Global"<<CenterPointRollGlobal<<std::endl; 00560 00561 LocalPoint CenterRollinDTFrame = DTSurfaceMB3.toLocal(CenterPointRollGlobal); 00562 //std::cout<<"\t \t \t Center (0,0,0) Roll In DT MB3 Local"<<CenterRollinDTFrame<<std::endl; 00563 00564 float D=CenterRollinDTFrame.z(); 00565 //std::cout<<"\t \t \t D="<<D<<"cm"<<std::endl; 00566 00567 X=Xo3+dx3*D/dz3; 00568 Y=Yo3+dy3*D/dz3; 00569 Z=D; 00570 00571 const RectangularStripTopology* top_=dynamic_cast<const RectangularStripTopology*>(&(rollasociated->topology())); 00572 LocalPoint xmin = top_->localPosition(0.); 00573 LocalPoint xmax = top_->localPosition((float)rollasociated->nstrips()); 00574 float rsize = fabs( xmax.x()-xmin.x() )*0.5; 00575 float stripl = top_->stripLength(); 00576 00577 //std::cout<<"\t \t \t X Predicted in DT MB3 Local= "<<X<<"cm"<<std::endl; 00578 //std::cout<<"\t \t \t Y Predicted in DT MB3 Local= "<<Y<<"cm"<<std::endl; 00579 //std::cout<<"\t \t \t Z Predicted in DT MB3 Local= "<<Z<<"cm"<<std::endl; 00580 00581 GlobalPoint GlobalPointExtrapolated = DTSurfaceMB3.toGlobal(LocalPoint(X,Y,Z)); 00582 //std::cout<<"\t \t \t Point ExtraPolated in Global"<<GlobalPointExtrapolated<< std::endl; 00583 00584 LocalPoint PointExtrapolatedRPCFrame = RPCSurfaceRB4.toLocal(GlobalPointExtrapolated); 00585 //std::cout<<"\t \t \t Point Extrapolated in RPC RB4 Local"<<PointExtrapolatedRPCFrame<< std::endl; 00586 00587 //std::cout<<"\t \t \t Does the extrapolation go inside this roll?"<<std::endl; 00588 //conditions to find the right roll to extrapolate 00589 if(fabs(PointExtrapolatedRPCFrame.z()) < 0.01 && 00590 fabs(PointExtrapolatedRPCFrame.x()) < rsize && 00591 fabs(PointExtrapolatedRPCFrame.y()) < stripl*0.5){ 00592 00593 //std::cout<<"\t \t \t \t yes"<<std::endl; 00594 //getting the number of the strip 00595 const float stripPredicted= 00596 rollasociated->strip(LocalPoint(PointExtrapolatedRPCFrame.x(),PointExtrapolatedRPCFrame.y(),0.)); 00597 00598 //std::cout<<"\t \t \t \t Candidate"<<rollasociated->id()<<" "<<"(from DT Segment) STRIP---> "<<stripPredicted<< std::endl; 00599 00600 //--------- HISTOGRAM STRIP PREDICTED FROM DT ------------------- 00601 00602 RPCDetId rollId = rollasociated->id(); 00603 uint32_t id = rollId.rawId(); 00604 00605 RPCDetId otherRollId1,otherRollId2; 00606 00607 if(rollId.roll() == 1){ 00608 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),2); 00609 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),3); 00610 otherRollId1 = tempRollId1; 00611 otherRollId2 = tempRollId2; 00612 } 00613 else if(rollId.roll() == 2){ 00614 00615 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),1); 00616 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),3); 00617 otherRollId1 = tempRollId1; 00618 otherRollId2 = tempRollId2; 00619 } 00620 else if(rollId.roll() == 3){ 00621 RPCDetId tempRollId1(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),1); 00622 RPCDetId tempRollId2(rollId.region(),rollId.ring(),rollId.station(),rollId.sector(),rollId.layer(),rollId.subsector(),2); 00623 otherRollId1 = tempRollId1; 00624 otherRollId2 = tempRollId2; 00625 } 00626 RPCDigiCollection::Range rpcRangeDigi1=rpcDigis->get(otherRollId1); 00627 RPCDigiCollection::Range rpcRangeDigi2=rpcDigis->get(otherRollId2); 00628 00629 _idList.push_back(id); 00630 00631 char detUnitLabel[128]; 00632 sprintf(detUnitLabel ,"%d",id); 00633 sprintf(layerLabel ,"layer%d_subsector%d_roll%d",rollId.layer(),rollId.subsector(),rollId.roll()); 00634 00635 std::map<uint32_t, std::map<std::string,MonitorElement*> >::iterator meItr = meCollection.find(id); 00636 if (meItr == meCollection.end()){ 00637 meCollection[id] = bookDetUnitMEEff(rollId); 00638 std::cout << "\t \t \t \t Create new histograms for "<<layerLabel<<std::endl; 00639 } 00640 00641 std::map<std::string, MonitorElement*> meMap=meCollection[id]; 00642 sprintf(meIdDT,"ExpectedOccupancyFromDT_%s",detUnitLabel); 00643 meMap[meIdDT]->Fill(stripPredicted); 00644 00645 sprintf(meIdDT,"ExpectedOccupancy2DFromDT_%s",detUnitLabel); 00646 meMap[meIdDT]->Fill(stripPredicted,Y); 00647 00648 //std::cout << "\t \t \t \t One for counterPREDICT"<<std::endl; 00649 totalcounter[0]++; 00650 buff=counter[0]; 00651 buff[rollId]++; 00652 counter[0]=buff; 00653 00654 bool anycoincidence=false; 00655 int stripDetected = 0; 00656 RPCDigiCollection::Range rpcRangeDigi = 00657 rpcDigis->get(rollasociated->id()); 00658 00659 int stripCounter = 0; 00660 00661 //loop over the digis in the event 00662 for (RPCDigiCollection::const_iterator digiIt = rpcRangeDigi.first;digiIt!=rpcRangeDigi.second;++digiIt){ 00663 stripCounter++; 00664 //std::cout<<"\t \t \t \t \t Digi "<<*digiIt<<std::endl;//print the digis in the event 00665 stripDetected=digiIt->strip(); 00666 00667 stripCounter++; 00668 float res = (float)(stripDetected) - stripPredicted; 00669 sprintf(meIdRPC,"RPCResiduals_%s",detUnitLabel); 00670 meMap[meIdRPC]->Fill(res); 00671 00672 sprintf(meIdRPC,"RPCResiduals2D_%s",detUnitLabel); 00673 meMap[meIdRPC]->Fill(res,Y); 00674 00675 if(res > 7) 00676 std::cout<<"STRANGE "<<"EVENTO NUM = "<<iEvent.id().event()<<" Residuo = "<<res<<" Strip Num = "<<stripDetected<<" Strip totali = "<<stripCounter<<std::endl; 00677 if(fabs((float)(stripDetected) - stripPredicted)<widestripsRB4){//Detected Vs Predicted 00678 //std::cout <<"\t \t \t \t \t COINCEDENCE Predict "<<stripPredicted<<" Detect "<<stripDetected<<std::endl; 00679 anycoincidence=true; 00680 //break;//funciona solo para hacerlo mas rapido 00681 } 00682 } 00683 if (anycoincidence){ 00684 00685 sprintf(meIdRPC,"ExpectedOccupancyFromDT_forCrT_%s",detUnitLabel); 00686 meMap[meIdRPC]->Fill(stripPredicted); 00687 sprintf(meIdRPC,"RealDetectedOccupancy_%s",detUnitLabel); 00688 meMap[meIdRPC]->Fill(stripDetected); 00689 sprintf(meIdRPC,"YExpectedOccupancyFromDT_%s",detUnitLabel); 00690 meMap[meIdRPC]->Fill(Y); 00691 00692 for (RPCDigiCollection::const_iterator digiIt1 = rpcRangeDigi1.first;digiIt1!=rpcRangeDigi1.second;++digiIt1){ 00693 00694 if(fabs(stripDetected - digiIt1->strip()) <= 1){ 00695 00696 sprintf(meIdRPC,"XCrossTalk_1_%s",detUnitLabel); 00697 meMap[meIdRPC]->Fill(stripPredicted); 00698 00699 sprintf(meIdRPC,"XDetectCrossTalk_1_%s",detUnitLabel); 00700 meMap[meIdRPC]->Fill(stripDetected); 00701 00702 sprintf(meIdRPC,"YCrossTalk_1_%s",detUnitLabel); 00703 meMap[meIdRPC]->Fill(Y); 00704 break; 00705 } 00706 } 00707 00708 for (RPCDigiCollection::const_iterator digiIt2 = rpcRangeDigi2.first;digiIt2!=rpcRangeDigi2.second;++digiIt2){ 00709 00710 if(fabs(stripDetected - digiIt2->strip()) <= 1){ 00711 00712 sprintf(meIdRPC,"XCrossTalk_2_%s",detUnitLabel); 00713 meMap[meIdRPC]->Fill(stripPredicted); 00714 00715 sprintf(meIdRPC,"XDetectCrossTalk_2_%s",detUnitLabel); 00716 meMap[meIdRPC]->Fill(stripDetected); 00717 00718 sprintf(meIdRPC,"YCrossTalk_2_%s",detUnitLabel); 00719 meMap[meIdRPC]->Fill(Y); 00720 break; 00721 } 00722 } 00723 00724 sprintf(meIdRPC,"RPCDataOccupancy_%s",detUnitLabel); 00725 meMap[meIdRPC]->Fill(stripPredicted); 00726 00727 sprintf(meIdRPC,"RPCDataOccupancy2D_%s",detUnitLabel); 00728 meMap[meIdRPC]->Fill(stripPredicted,Y); 00729 00730 totalcounter[1]++; 00731 buff=counter[1]; 00732 buff[rollId]++; 00733 counter[1]=buff; 00734 } 00735 else{ 00736 totalcounter[2]++; 00737 buff=counter[2]; 00738 buff[rollId]++; 00739 counter[2]=buff; 00740 //std::cout <<"\t \t \t \t \t XXXXX THIS PREDICTION DOESN'T HAVE ANY CORRESPONDENCE WITH THE DATA"<<std::endl; 00741 //std::cout << "\t \t \t \t \t One for counterFAIL"<<std::endl; 00742 //ofrej<<"Wh "<<dtWheel<<"\t St "<<dtStation 00743 // <<"\t Se "<<dtSector<<"\t Event " 00744 // <<iEvent.id().event()<<std::endl; 00745 00746 } 00747 } 00748 else{ 00749 //std::cout<<"\t \t \t \t no"<<std::endl; 00750 }//Condition for the right match 00751 00752 }//loop over all the rolls FOR RB3 ------------------------------- 00753 }//Si tenemos solo un segmento en la chamber 00754 }//Condition Segment in MB3 00755 }//Loop avoer all the segments to see if it is the right MB3 00756 } 00757 else{ 00758 //std::cout<<"\t \t Strange Segment"<<std::endl; 00759 }//Is not a 4D Segment neither a 2D in MB4 00760 } 00761 else { 00762 //std::cout<<"\t \t no"<<std::endl; 00763 }//There is one segment in the chamber? 00764 }//loop over the segments 00765 }else { 00766 //std::cout<<"This Event doesn't have any DT4DSegment"<<std::endl; 00767 }//is ther more than 1 segment in this event? 00768 00769 }
std::map< std::string, MonitorElement * > RPCMonitorEfficiency::bookDetUnitMEEff | ( | RPCDetId & | detId | ) |
Booking of MonitoringElemnt for one RPCDetId (= roll).
Name components common to current RPDDetId
Definition at line 15 of file RPCBookDetUnitMEEff.cc.
References DQMStore::book1D(), DQMStore::book2D(), dbe, detId, RPCDetId::layer(), RPCDetId::region(), RPCDetId::ring(), RPCDetId::roll(), RPCDetId::sector(), DQMStore::setCurrentFolder(), RPCDetId::station(), and RPCDetId::subsector().
Referenced by analyze().
00015 { 00016 00017 std::map<std::string, MonitorElement*> meMap; 00018 std::string regionName; 00019 std::string ringType; 00020 if(detId.region()==0) { 00021 regionName="Barrel"; 00022 ringType="Wheel"; 00023 } 00024 else{ 00025 ringType="Disk"; 00026 if(detId.region() == -1) regionName="Encap-"; 00027 if(detId.region() == 1) regionName="Encap+"; 00028 } 00029 00030 char folder[120]; 00031 sprintf(folder,"RPC/Efficiency/%s/%s_%d/station_%d/sector_%d",regionName.c_str(),ringType.c_str(), 00032 detId.ring(),detId.station(),detId.sector()); 00033 00034 //std::cout<<"BOOKING 1"<<std::endl; 00035 dbe->setCurrentFolder(folder); 00036 //std::cout<<"BOOKING 2"<<std::endl; 00038 char detUnitLabel[128]; 00039 char layerLabel[128]; 00040 sprintf(detUnitLabel ,"%d",detId()); 00041 sprintf(layerLabel ,"layer%d_subsector%d_roll%d",detId.layer(),detId.subsector(),detId.roll()); 00042 //std::cout<<"BOOKING 3"<<std::endl; 00043 00044 char meId [128]; 00045 char meTitle [128]; 00046 00047 // Begin booking 00048 sprintf(meId,"ExpectedOccupancyFromDT_%s",detUnitLabel); 00049 sprintf(meTitle,"ExpectedOccupancyFromDT_for_%s",layerLabel); 00050 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00051 00052 sprintf(meId,"ExpectedOccupancyFromDT_forCrT_%s",detUnitLabel); 00053 sprintf(meTitle,"ExpectedOccupancyFromDT_forCrT_for_%s",layerLabel); 00054 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00055 00056 sprintf(meId,"YExpectedOccupancyFromDT_%s",detUnitLabel); 00057 sprintf(meTitle,"YExpectedOccupancyFromDT_for_%s",layerLabel); 00058 meMap[meId] = dbe->book1D(meId, meTitle, 100, -100, 100); 00059 00060 sprintf(meId,"ExpectedOccupancy2DFromDT_%s",detUnitLabel); 00061 sprintf(meTitle,"ExpectedOccupancy2DFromDT_for_%s",layerLabel); 00062 meMap[meId] = dbe->book2D(meId, meTitle, 100, 0.5, 100.5,200,-100.,100.); 00063 00064 //std::cout<<"BOOKING 4"<<std::endl; 00065 sprintf(meId,"RPCDataOccupancy_%s",detUnitLabel); 00066 sprintf(meTitle,"RPCDataOccupancy_for_%s",layerLabel); 00067 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00068 00069 sprintf(meId,"RPCDataOccupancy2D_%s",detUnitLabel); 00070 sprintf(meTitle,"RPCDataOccupancy2D_for_%s",layerLabel); 00071 meMap[meId] = dbe->book2D(meId, meTitle, 100, 0.5, 100.5,200,-100.,100.); 00072 00073 sprintf(meId,"RealDetectedOccupancy_%s",detUnitLabel); 00074 sprintf(meTitle,"RealDetectedOccupancy_for_%s",layerLabel); 00075 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00076 00077 //std::cout<<"BOOKING 5"<<std::endl; 00078 sprintf(meId,"EfficienyFromDTExtrapolation_%s",detUnitLabel); 00079 sprintf(meTitle,"EfficienyFromDTExtrapolation_for_%s",layerLabel); 00080 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00081 00082 sprintf(meId,"EfficienyFromDT2DExtrapolation_%s",detUnitLabel); 00083 sprintf(meTitle,"EfficienyFromDT2DExtrapolation_for_%s",layerLabel); 00084 meMap[meId] = dbe->book2D(meId, meTitle, 100, 0.5, 100.5,200,-100.,100.); 00085 00086 sprintf(meId,"XCrossTalkFromDTExtrapolation_1_%s",detUnitLabel); 00087 sprintf(meTitle,"XCrossTalkFromDTExtrapolation_for_%s",layerLabel); 00088 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00089 00090 sprintf(meId,"XCrossTalkFromDTExtrapolation_2_%s",detUnitLabel); 00091 sprintf(meTitle,"XCrossTalkFromDTExtrapolation_for_%s",layerLabel); 00092 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00093 00094 sprintf(meId,"XCrossTalkFromDetectedStrip_1_%s",detUnitLabel); 00095 sprintf(meTitle,"XCrossTalkFromDetectedStrip_1_for_%s",layerLabel); 00096 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00097 00098 sprintf(meId,"XCrossTalkFromDetectedStrip_2_%s",detUnitLabel); 00099 sprintf(meTitle,"XCrossTalkFromDetectedStrip_2_for_%s",layerLabel); 00100 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00101 00102 sprintf(meId,"YCrossTalkFromDTExtrapolation_1_%s",detUnitLabel); 00103 sprintf(meTitle,"YCrossTalkFromDTExtrapolation_1_for_%s",layerLabel); 00104 meMap[meId] = dbe->book1D(meId, meTitle, 200, -100., 100.); 00105 00106 sprintf(meId,"YCrossTalkFromDTExtrapolation_2_%s",detUnitLabel); 00107 sprintf(meTitle,"YCrossTalkFromDTExtrapolation_2_for_%s",layerLabel); 00108 meMap[meId] = dbe->book1D(meId, meTitle, 200, -100., 100.); 00109 00110 sprintf(meId,"RPCResiduals_%s",detUnitLabel); 00111 sprintf(meTitle,"RPCResiduals_for_%s",layerLabel); 00112 meMap[meId] = dbe->book1D(meId, meTitle, 201,-100.5, 100.5); 00113 00114 sprintf(meId,"RPCResiduals2D_%s",detUnitLabel); 00115 sprintf(meTitle,"RPCResiduals2D_for_%s",layerLabel); 00116 meMap[meId] = dbe->book2D(meId, meTitle, 201,-100.5, 100.5,200,-100.,100.); 00117 00118 sprintf(meId,"XCrossTalk_1_%s",detUnitLabel); 00119 sprintf(meTitle,"XCrossTalk_1_for_%s",layerLabel); 00120 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00121 00122 sprintf(meId,"XCrossTalk_2_%s",detUnitLabel); 00123 sprintf(meTitle,"XCrossTalk_2_for_%s",layerLabel); 00124 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00125 00126 sprintf(meId,"XDetectCrossTalk_1_%s",detUnitLabel); 00127 sprintf(meTitle,"XDetectCrossTalk_1_for_%s",layerLabel); 00128 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00129 00130 sprintf(meId,"XDetectCrossTalk_2_%s",detUnitLabel); 00131 sprintf(meTitle,"XDetectCrossTalk_2_for_%s",layerLabel); 00132 meMap[meId] = dbe->book1D(meId, meTitle, 100, 0.5, 100.5); 00133 00134 sprintf(meId,"YCrossTalk_1_%s",detUnitLabel); 00135 sprintf(meTitle,"YCrossTalk_1_for_%s",layerLabel); 00136 meMap[meId] = dbe->book1D(meId, meTitle, 100, -100, 100); 00137 00138 sprintf(meId,"YCrossTalk_2_%s",detUnitLabel); 00139 sprintf(meTitle,"YCrossTalk_2_for_%s",layerLabel); 00140 meMap[meId] = dbe->book1D(meId, meTitle, 200, -100., 100.); 00141 00142 return meMap; 00143 }
Reimplemented from edm::EDAnalyzer.
Definition at line 771 of file RPCMonitorEfficiency.cc.
References _idList, counter, GenMuonPlsPt100GeV_cfg::cout, dbe, EffRootFileName, EffSaveRootFile, lat::endl(), i, j, meCollection, p, r, DQMStore::save(), funct::sqrt(), and totalcounter.
00771 { 00772 std::map<RPCDetId, int> pred = counter[0]; 00773 std::map<RPCDetId, int> obse = counter[1]; 00774 std::map<RPCDetId, int> reje = counter[2]; 00775 00776 std::map<RPCDetId, int>::iterator irpc; 00777 for (irpc=pred.begin(); irpc!=pred.end();irpc++){ 00778 RPCDetId id=irpc->first; 00779 int p=pred[id]; 00780 int o=obse[id]; 00781 int r=reje[id]; 00782 assert(p==o+r); 00783 float ef = float(o)/float(p); 00784 float er = sqrt(ef*(1.-ef)/float(p)); 00785 std::cout <<"\n "<<id<<"\t Predicted "<<p<<"\t Observed "<<o<<"\t Eff = "<<ef*100.<<" % +/- "<<er*100.<<" %"; 00786 if(ef<0.8){ 00787 std::cout<<"\t \t Warning!"; 00788 } 00789 } 00790 00791 float tote = float(totalcounter[1])/float(totalcounter[0]); 00792 float totr = sqrt(tote*(1.-tote)/float(totalcounter[0])); 00793 std::cout <<"\n\n \t \t TOTAL EFFICIENCY \t Predicted "<<totalcounter[1]<<"\t Observed "<<totalcounter[0]<<"\t Eff = "<<tote*100.<<"\t +/- \t"<<totr*100.<<"%"<<std::endl; 00794 std::cout <<totalcounter[1]<<" "<<totalcounter[0]<<" flagcode"<<std::endl; 00795 00796 00797 std::vector<uint32_t>::iterator meIt; 00798 for(meIt = _idList.begin(); meIt != _idList.end(); ++meIt){ 00799 00800 char detUnitLabel[128]; 00801 char meIdRPC [128]; 00802 char meIdDT [128]; 00803 char meIdDTCrT [128]; 00804 char effIdRPC [128]; 00805 char meIdRPCCrT [128]; 00806 00807 char XCrTalkId_1 [128]; 00808 char XCrTalkId_2 [128]; 00809 char effXCrTalkId_1 [128]; 00810 char effXCrTalkId_2 [128]; 00811 00812 char XCrTalkDetId_1 [128]; 00813 char XCrTalkDetId_2 [128]; 00814 char effXCrTalkDetId_1 [128]; 00815 char effXCrTalkDetId_2 [128]; 00816 00817 char YCrTalkId_1 [128]; 00818 char YCrTalkId_2 [128]; 00819 char effYCrTalkId_1 [128]; 00820 char effYCrTalkId_2 [128]; 00821 00822 char YCrTalkDTId [128]; 00823 00824 char meIdRPC_2D [128]; 00825 char meIdDT_2D [128]; 00826 char effIdRPC_2D [128]; 00827 00828 sprintf(detUnitLabel ,"%d",*meIt); 00829 sprintf(meIdRPC,"RPCDataOccupancy_%s",detUnitLabel); 00830 sprintf(meIdDT,"ExpectedOccupancyFromDT_%s",detUnitLabel); 00831 sprintf(meIdDTCrT,"ExpectedOccupancyFromDT_forCrT_%s",detUnitLabel); 00832 sprintf(YCrTalkDTId,"YExpectedOccupancyFromDT_%s",detUnitLabel); 00833 sprintf(meIdRPCCrT,"RealDetectedOccupancy_%s",detUnitLabel); 00834 00835 sprintf(effIdRPC,"EfficienyFromDTExtrapolation_%s",detUnitLabel); 00836 00837 sprintf(meIdRPC_2D,"RPCDataOccupancy2D_%s",detUnitLabel); 00838 sprintf(meIdDT_2D,"ExpectedOccupancy2DFromDT_%s",detUnitLabel); 00839 sprintf(effIdRPC_2D,"EfficienyFromDT2DExtrapolation_%s",detUnitLabel); 00840 00841 sprintf(XCrTalkId_1,"XCrossTalk_1_%s",detUnitLabel); 00842 sprintf(XCrTalkId_2,"XCrossTalk_2_%s",detUnitLabel); 00843 sprintf(YCrTalkId_1,"YCrossTalk_1_%s",detUnitLabel); 00844 sprintf(YCrTalkId_2,"YCrossTalk_2_%s",detUnitLabel); 00845 sprintf(XCrTalkDetId_1,"XDetectCrossTalk_1_%s",detUnitLabel); 00846 sprintf(XCrTalkDetId_2,"XDetectCrossTalk_2_%s",detUnitLabel); 00847 00848 sprintf(effXCrTalkId_1,"XCrossTalkFromDTExtrapolation_1_%s",detUnitLabel); 00849 sprintf(effXCrTalkId_2,"XCrossTalkFromDTExtrapolation_2_%s",detUnitLabel); 00850 sprintf(effXCrTalkDetId_1,"XCrossTalkFromDetectedStrip_1_%s",detUnitLabel); 00851 sprintf(effXCrTalkDetId_2,"XCrossTalkFromDetectedStrip_2_%s",detUnitLabel); 00852 sprintf(effYCrTalkId_1,"YCrossTalkFromDTExtrapolation_1_%s",detUnitLabel); 00853 sprintf(effYCrTalkId_2,"YCrossTalkFromDTExtrapolation_2_%s",detUnitLabel); 00854 00855 std::map<std::string, MonitorElement*> meMap=meCollection[*meIt]; 00856 00857 for(unsigned int i=1;i<=100;++i){ 00858 00859 if(meMap[meIdDT]->getBinContent(i) != 0){ 00860 float eff = meMap[meIdRPC]->getBinContent(i)/meMap[meIdDT]->getBinContent(i); 00861 float erreff = sqrt(eff*(1-eff)/meMap[meIdDT]->getBinContent(i)); 00862 meMap[effIdRPC]->setBinContent(i,eff*100.); 00863 meMap[effIdRPC]->setBinError(i,erreff*100.); 00864 } 00865 } 00866 for(unsigned int i=1;i<=100;++i){ 00867 for(unsigned int j=1;j<=200;++j){ 00868 if(meMap[meIdDT_2D]->getBinContent(i,j) != 0){ 00869 float eff = meMap[meIdRPC_2D]->getBinContent(i,j)/meMap[meIdDT_2D]->getBinContent(i,j); 00870 float erreff = sqrt(eff*(1-eff)/meMap[meIdDT_2D]->getBinContent(i,j)); 00871 meMap[effIdRPC_2D]->setBinContent(i,j,eff*100.); 00872 meMap[effIdRPC_2D]->setBinError(i,j,erreff*100.); 00873 } 00874 } 00875 } 00876 00877 //-------------------- CROSS TALK PLOT --------------------------------------------------- 00878 00879 //-------------------- With predicted STRIP ----------------------------------------------- 00880 00881 for(unsigned int i=1;i<=100;++i){ 00882 00883 if(meMap[meIdDTCrT]->getBinContent(i) != 0){ 00884 float crt = meMap[XCrTalkId_1]->getBinContent(i)/meMap[meIdDTCrT]->getBinContent(i); 00885 float errcrt = sqrt(crt*(1-crt)/meMap[meIdDTCrT]->getBinContent(i)); 00886 meMap[effXCrTalkId_1]->setBinContent(i,crt*100.); 00887 meMap[effXCrTalkId_1]->setBinError(i,errcrt*100.); 00888 } 00889 } 00890 00891 for(unsigned int i=1;i<=100;++i){ 00892 00893 if(meMap[meIdDTCrT]->getBinContent(i) != 0){ 00894 float crt = meMap[XCrTalkId_2]->getBinContent(i)/meMap[meIdDTCrT]->getBinContent(i); 00895 float errcrt = sqrt(crt*(1-crt)/meMap[meIdDTCrT]->getBinContent(i)); 00896 meMap[effXCrTalkId_2]->setBinContent(i,crt*100.); 00897 meMap[effXCrTalkId_2]->setBinError(i,errcrt*100.); 00898 } 00899 } 00900 00901 //-------------------- With detected STRIP ----------------------------------------------- 00902 00903 for(unsigned int i=1;i<=100;++i){ 00904 if(meMap[meIdRPCCrT]->getBinContent(i) != 0){ 00905 float crt = meMap[XCrTalkDetId_1]->getBinContent(i)/meMap[meIdRPCCrT]->getBinContent(i); 00906 float errcrt = sqrt(crt*(1-crt)/meMap[meIdRPCCrT]->getBinContent(i)); 00907 meMap[effXCrTalkDetId_1]->setBinContent(i,crt*100.); 00908 meMap[effXCrTalkDetId_1]->setBinError(i,errcrt*100.); 00909 } 00910 } 00911 00912 for(unsigned int i=1;i<=100;++i){ 00913 00914 if(meMap[meIdRPCCrT]->getBinContent(i) != 0){ 00915 float crt = meMap[XCrTalkDetId_2]->getBinContent(i)/meMap[meIdRPCCrT]->getBinContent(i); 00916 float errcrt = sqrt(crt*(1-crt)/meMap[meIdRPCCrT]->getBinContent(i)); 00917 meMap[effXCrTalkDetId_2]->setBinContent(i,crt*100.); 00918 meMap[effXCrTalkDetId_2]->setBinError(i,errcrt*100.); 00919 } 00920 } 00921 00922 //-------------------- With Y coordinate ------------------------------------------------- 00923 00924 for(unsigned int i=1;i<=200;++i){ 00925 00926 if(meMap[YCrTalkDTId]->getBinContent(i) != 0){ 00927 float crt = meMap[YCrTalkId_1]->getBinContent(i)/meMap[YCrTalkDTId]->getBinContent(i); 00928 float errcrt = sqrt(crt*(1-crt)/meMap[YCrTalkDTId]->getBinContent(i)); 00929 meMap[effYCrTalkId_1]->setBinContent(i,crt*100.); 00930 meMap[effYCrTalkId_1]->setBinError(i,errcrt*100.); 00931 } 00932 } 00933 00934 for(unsigned int i=1;i<=200;++i){ 00935 00936 if(meMap[YCrTalkDTId]->getBinContent(i) != 0){ 00937 float crt = meMap[YCrTalkId_2]->getBinContent(i)/meMap[YCrTalkDTId]->getBinContent(i); 00938 float errcrt = sqrt(crt*(1-crt)/meMap[YCrTalkDTId]->getBinContent(i)); 00939 meMap[effYCrTalkId_2]->setBinContent(i,crt*100.); 00940 meMap[effYCrTalkId_2]->setBinError(i,errcrt*100.); 00941 } 00942 } 00943 } 00944 00945 if(EffSaveRootFile) dbe->save(EffRootFileName); 00946 // theFile->Write(); 00947 // theile->Close(); 00948 std::cout<<"End Job"<<std::endl; 00949 }
std::vector<uint32_t> RPCMonitorEfficiency::_idList [private] |
Definition at line 68 of file RPCMonitorEfficiency.h.
Referenced by analyze(), endJob(), and RPCMonitorEfficiency().
std::vector<std::map<RPCDetId, int> > RPCMonitorEfficiency::counter [private] |
Definition at line 70 of file RPCMonitorEfficiency.h.
Referenced by analyze(), endJob(), and RPCMonitorEfficiency().
DQMStore* RPCMonitorEfficiency::dbe [private] |
back-end interface
Definition at line 59 of file RPCMonitorEfficiency.h.
Referenced by bookDetUnitMEEff(), endJob(), and RPCMonitorEfficiency().
bool RPCMonitorEfficiency::debug [private] |
Definition at line 62 of file RPCMonitorEfficiency.h.
std::string RPCMonitorEfficiency::digiLabel [private] |
Definition at line 64 of file RPCMonitorEfficiency.h.
Referenced by analyze(), and RPCMonitorEfficiency().
std::string RPCMonitorEfficiency::EffRootFileName [private] |
Definition at line 55 of file RPCMonitorEfficiency.h.
Referenced by endJob(), and RPCMonitorEfficiency().
bool RPCMonitorEfficiency::EffSaveRootFile [private] |
Definition at line 53 of file RPCMonitorEfficiency.h.
Referenced by endJob(), and RPCMonitorEfficiency().
MonitorElement* RPCMonitorEfficiency::h1 [private] |
Definition at line 60 of file RPCMonitorEfficiency.h.
TH1F* RPCMonitorEfficiency::hPositionX [private] |
Definition at line 67 of file RPCMonitorEfficiency.h.
std::map<uint32_t, std::map<std::string, MonitorElement*> > RPCMonitorEfficiency::meCollection [private] |
std::string RPCMonitorEfficiency::nameInLog [private] |
Definition at line 52 of file RPCMonitorEfficiency.h.
std::ofstream RPCMonitorEfficiency::ofrej [private] |
Definition at line 72 of file RPCMonitorEfficiency.h.
std::string RPCMonitorEfficiency::theRecHits4DLabel [private] |
Definition at line 63 of file RPCMonitorEfficiency.h.
Referenced by analyze(), and RPCMonitorEfficiency().
std::vector<int> RPCMonitorEfficiency::totalcounter [private] |
Definition at line 71 of file RPCMonitorEfficiency.h.
Referenced by analyze(), endJob(), and RPCMonitorEfficiency().