Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:27:34

0001 
0002 #include "TVector3.h"
0003 #include "TH2.h"
0004 #include "TLine.h"
0005 #include "TLatex.h"
0006 #include "TPaveText.h"
0007 #include "TCanvas.h"
0008 #include "TEveWindow.h"
0009 
0010 #include "DataFormats/TrackReco/interface/Track.h"
0011 #include "DataFormats/TrackReco/interface/HitPattern.h"
0012 #include "DataFormats/TrackingRecHit/interface/TrackingRecHit.h"
0013 
0014 #include "Fireworks/Core/interface/FWDetailView.h"
0015 #include "Fireworks/Core/interface/FWGeometry.h"
0016 #include "Fireworks/Core/interface/FWEventItem.h"
0017 #include "Fireworks/Core/interface/FWModelId.h"
0018 #include "Fireworks/Core/interface/fwLog.h"
0019 #include "Fireworks/Tracks/plugins/FWTrackResidualDetailView.h"
0020 
0021 using reco::HitPattern;
0022 using reco::Track;
0023 using reco::TrackBase;
0024 using reco::TrackResiduals;
0025 
0026 const char* FWTrackResidualDetailView::m_det_tracker_str[6] = {"PXB", "PXF", "TIB", "TID", "TOB", "TEC"};
0027 
0028 FWTrackResidualDetailView::FWTrackResidualDetailView()
0029     : m_ndet(0),
0030       m_nhits(0),
0031       m_resXFill(3007),
0032       m_resXCol(kGreen - 9),
0033       m_resYFill(3006),
0034       m_resYCol(kWhite),
0035       m_stereoFill(3004),
0036       m_stereoCol(kCyan - 9),
0037       m_invalidFill(3001),
0038       m_invalidCol(kRed) {
0039   memset(m_det, 0, sizeof(m_det));
0040   memset(res, 0, sizeof(res));
0041   memset(hittype, 0, sizeof(hittype));
0042   memset(stereo, 0, sizeof(stereo));
0043   memset(substruct, 0, sizeof(substruct));
0044   memset(subsubstruct, 0, sizeof(subsubstruct));
0045   memset(m_detector, 0, sizeof(m_detector));
0046 }
0047 
0048 FWTrackResidualDetailView::~FWTrackResidualDetailView() {}
0049 
0050 void FWTrackResidualDetailView::prepareData(const FWModelId& id, const reco::Track* track) {
0051   auto const& residuals = track->residuals();
0052 
0053   const FWGeometry* geom = id.item()->getGeom();
0054   assert(geom != nullptr);
0055 
0056   const HitPattern& hitpat = track->hitPattern();
0057   m_nhits = hitpat.numberOfAllHits(reco::HitPattern::TRACK_HITS);
0058   for (int i = 0; i < m_nhits; ++i) {
0059     //printf("there are %d hits in the pattern, %d in the vector, this is %u\n",
0060     //        m_nhits, track->recHitsEnd() - track->recHitsBegin(), (*(track->recHitsBegin() + i))->geographicalId().rawId());
0061     uint32_t pattern = hitpat.getHitPattern(reco::HitPattern::TRACK_HITS, i);
0062     hittype[i] = HitPattern::getHitType(pattern);
0063     stereo[i] = HitPattern::getSide(pattern);
0064     subsubstruct[i] = HitPattern::getSubSubStructure(pattern);
0065     substruct[i] = HitPattern::getSubStructure(pattern);
0066     m_detector[i] = HitPattern::getSubDetector(pattern);
0067     if ((*(track->recHitsBegin() + i))->isValid()) {
0068       res[0][i] = getSignedResidual(geom, (*(track->recHitsBegin() + i))->geographicalId().rawId(), residuals.pullX(i));
0069     } else {
0070       res[0][i] = 0;
0071     }
0072     res[1][i] = residuals.pullY(i);
0073     // printf("%s, %i\n",m_det_tracker_str[substruct[i]-1],subsubstruct[i]);
0074   }
0075 
0076   m_det[0] = 0;
0077   for (int j = 0; j < m_nhits - 1;) {
0078     int k = j + 1;
0079     for (; k < m_nhits; k++) {
0080       if (substruct[j] == substruct[k] && subsubstruct[j] == subsubstruct[k]) {
0081         if (k == (m_nhits - 1))
0082           j = k;
0083       } else {
0084         m_ndet++;
0085         j = k;
0086         m_det[m_ndet] = j;
0087         break;
0088       }
0089     }
0090   }
0091   m_ndet++;
0092   m_det[m_ndet] = m_nhits;
0093   // printDebug();
0094 }
0095 
0096 void FWTrackResidualDetailView::build(const FWModelId& id, const reco::Track* track) {
0097   if (!track->extra().isAvailable()) {
0098     fwLog(fwlog::kError) << " no track extra information is available.\n";
0099     m_viewCanvas->cd();
0100     TLatex* latex = new TLatex();
0101     latex->SetTextSize(0.1);
0102     latex->DrawLatex(0.1, 0.5, "No track extra information");
0103     return;
0104   }
0105   prepareData(id, track);
0106 
0107   // draw histogram
0108   m_viewCanvas->cd();
0109   m_viewCanvas->SetHighLightColor(-1);
0110   TH2F* h_res = new TH2F("h_resx", "h_resx", 10, -5.5, 5.5, m_ndet, 0, m_ndet);
0111   TPad* padX = new TPad("pad1", "pad1", 0.2, 0., 0.8, 0.99);
0112   padX->SetBorderMode(0);
0113   padX->SetLeftMargin(0.2);
0114   padX->Draw();
0115   padX->cd();
0116   padX->SetFrameLineWidth(0);
0117   padX->Modified();
0118   h_res->SetDirectory(nullptr);
0119   h_res->SetStats(kFALSE);
0120   h_res->SetTitle("");
0121   h_res->SetXTitle("residual");
0122   h_res->GetXaxis()->SetTickLength(0);
0123   h_res->GetYaxis()->SetTickLength(0);
0124   h_res->GetXaxis()->SetNdivisions(20);
0125   h_res->GetYaxis()->SetLabelSize(0.06);
0126   h_res->Draw();
0127   padX->SetGridy();
0128 
0129   float larray[9] = {0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.5, 4.5, 5.5};
0130   float larray2[8];
0131   for (int l = 0; l < 8; l++) {
0132     float diff2 = (larray[l + 1] - larray[l]) / 2;
0133     larray2[l] = larray[l] + diff2;
0134     //  printf("(%.1f,%.1f),",larray[i],larray[i+1]);
0135   }
0136 
0137   int resi[2][64];
0138   for (int l = 0; l < m_nhits; l++) {
0139     for (int k = 0; k < 8; k++) {
0140       if (fabs(res[0][l]) == larray2[k])
0141         resi[0][l] = k;
0142       if (fabs(res[1][l]) == larray2[k])
0143         resi[1][l] = k;
0144     }
0145   }
0146 
0147   TLine* lines[17];
0148   for (int l = 0; l < 17; l++) {
0149     int ix = l % 9;
0150     int sign = 1;
0151     sign = (l > 8) ? -1 : 1;
0152     lines[l] = new TLine(sign * larray[ix], 0, sign * larray[ix], m_ndet);
0153     if (l != 9)
0154       lines[l]->SetLineStyle(3);
0155     padX->cd();
0156     lines[l]->Draw();
0157   }
0158 
0159   float width = 0.25;
0160   int filltype;
0161   Color_t color;
0162   Double_t box[4];
0163   padX->cd();
0164 
0165   for (int h = 0; h < 2; h++) {
0166     float height1 = 0;
0167     for (int j = 0; j < m_ndet; j++) {
0168       // take only X res and Y pixel residals
0169       if (strcmp(m_det_tracker_str[substruct[m_det[j]] - 1], "PXB") && h)
0170         continue;
0171 
0172       char det_str2[256];
0173       snprintf(det_str2, 255, "%s/%i", m_det_tracker_str[substruct[m_det[j]] - 1], subsubstruct[m_det[j]]);
0174       h_res->GetYaxis()->SetBinLabel(j + 1, det_str2);
0175 
0176       int diff = m_det[j + 1] - m_det[j];
0177       int k = 0;
0178       width = 1.0 / diff;
0179 
0180       for (int l = m_det[j]; l < (m_det[j] + diff); l++) {
0181         //      g->SetPoint(l,resx[l],j+0.5);
0182         //  printf("%i, %f %f %f\n",l,resx[l],sign*larray[resxi[l]],sign*larray[resxi[l]+1]);
0183         int sign = (res[h][l] < 0) ? -1 : 1;
0184         box[0] = (hittype[l] == 0) ? sign * larray[resi[h][l]] : -5.5;
0185         box[2] = (hittype[l] == 0) ? sign * larray[resi[h][l] + 1] : 5.5;
0186         box[1] = height1 + width * k;
0187         box[3] = height1 + width * (k + 1);
0188 
0189         if (stereo[l] == 1) {
0190           color = m_stereoCol;
0191           filltype = m_stereoFill;
0192         } else if (hittype[l] != 0) {
0193           color = m_invalidCol;
0194           filltype = m_invalidFill;
0195         } else {
0196           filltype = h ? m_resYFill : m_resXFill;
0197           color = h ? m_resYCol : m_resXCol;
0198         }
0199 
0200         drawCanvasBox(box, color, filltype, h < 1);
0201         k++;
0202       }
0203       height1 += 1;
0204     }
0205   }
0206 
0207   //  title
0208   const char* res_str = "residuals in Si detector local x-y coord.";
0209   TPaveText* pt = new TPaveText(0.0, 0.91, 1, 0.99, "blNDC");
0210   pt->SetBorderSize(0);
0211   pt->SetFillColor(kWhite);
0212   pt->AddText(res_str);
0213   pt->Draw();
0214 
0215   m_viewCanvas->cd();
0216   m_viewCanvas->SetEditable(kFALSE);
0217 
0218   setTextInfo(id, track);
0219 }
0220 
0221 double FWTrackResidualDetailView::getSignedResidual(const FWGeometry* geom, unsigned int id, double resX) {
0222   double local1[3] = {0, 0, 0};
0223   double local2[3] = {resX, 0, 0};
0224   double global1[3], global2[3];
0225   const TGeoMatrix* m = geom->getMatrix(id);
0226   assert(m != nullptr);
0227   m->LocalToMaster(local1, global1);
0228   m->LocalToMaster(local2, global2);
0229   TVector3 g1 = global1;
0230   TVector3 g2 = global2;
0231   if (g2.DeltaPhi(g1) > 0)
0232     return resX;
0233   else
0234     return -resX;
0235 }
0236 
0237 void FWTrackResidualDetailView::printDebug() {
0238   for (int i = 0; i < m_ndet; i++) {
0239     std::cout << "<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<" << std::endl;
0240     std::cout << "idx " << i << " det[idx] " << m_det[i] << std::endl;
0241     std::cout << "m_det idx " << m_det[i] << std::endl;
0242     std::cout << "m_det_tracker_str idx " << substruct[m_det[i]] - 1 << std::endl;
0243     printf("m_det[idx] %i m_det_tracker_str %s substruct %i\n",
0244            m_det[i],
0245            m_det_tracker_str[substruct[m_det[i]] - 1],
0246            subsubstruct[m_det[i]]);
0247   }
0248 }
0249 
0250 void FWTrackResidualDetailView::setTextInfo(const FWModelId& /*id*/, const reco::Track* /*track*/) {
0251   m_infoCanvas->cd();
0252 
0253   char mytext[256];
0254   Double_t fontsize = 0.07;
0255   TLatex* latex = new TLatex();
0256   latex->SetTextSize(fontsize);
0257   latex->Draw();
0258   Double_t x0 = 0.02;
0259   Double_t y = 0.95;
0260 
0261   // summary
0262   int nvalid = 0;
0263   int npix = 0;
0264   int nstrip = 0;
0265   for (int i = 0; i < m_nhits; i++) {
0266     if (hittype[i] == 0)
0267       nvalid++;
0268     if (substruct[i] < 3)
0269       npix++;
0270     else
0271       nstrip++;
0272   }
0273 
0274   latex->SetTextSize(fontsize);
0275   Double_t boxH = 0.25 * fontsize;
0276 
0277   double yStep = 0.04;
0278 
0279   latex->DrawLatex(x0, y, "Residual:");
0280   y -= yStep;
0281   latex->DrawLatex(
0282       x0,
0283       y,
0284       "sgn(#hat{X}#bullet#hat{#phi}) #times #frac{X_{hit} - X_{traj}}{#sqrt{#sigma^{2}_{hit} + #sigma^{2}_{traj}}}");
0285   y -= 2.5 * yStep;
0286   snprintf(mytext, 255, "layers hit: %i", m_ndet);
0287   latex->DrawLatex(x0, y, mytext);
0288   y -= yStep;
0289   snprintf(mytext, 255, "valid Si hits: %i", nvalid);
0290   latex->DrawLatex(x0, y, mytext);
0291   y -= yStep;
0292   snprintf(mytext, 255, "total Si hits: %i", m_nhits);
0293   latex->DrawLatex(x0, y, mytext);
0294   y -= yStep;
0295   snprintf(mytext, 255, "valid Si pixel hits: %i", npix);
0296   latex->DrawLatex(x0, y, mytext);
0297   y -= yStep;
0298   snprintf(mytext, 255, "valid Si strip hits: %i", nstrip);
0299   latex->DrawLatex(x0, y, mytext);
0300 
0301   Double_t pos[4];
0302   pos[0] = 0.4;
0303   pos[2] = 0.55;
0304 
0305   y -= yStep * 2;
0306   latex->DrawLatex(x0, y, "X hit");
0307   pos[1] = y;
0308   pos[3] = pos[1] + boxH;
0309   drawCanvasBox(pos, m_resXCol, m_resXFill);
0310 
0311   y -= yStep;
0312   latex->DrawLatex(x0, y, "Y hit");
0313   pos[1] = y;
0314   pos[3] = pos[1] + boxH;
0315   drawCanvasBox(pos, m_resYCol, m_resYFill, false);
0316 
0317   y -= yStep;
0318   latex->DrawLatex(x0, y, "stereo hit");
0319   pos[1] = y;
0320   pos[3] = pos[1] + boxH;
0321   drawCanvasBox(pos, m_stereoCol, m_stereoFill);
0322 
0323   y -= yStep;
0324   latex->DrawLatex(x0, y, "invalid hit");
0325   pos[1] = y;
0326   pos[3] = pos[1] + boxH;
0327   drawCanvasBox(pos, m_invalidCol, m_invalidFill);
0328 }
0329 
0330 REGISTER_FWDETAILVIEW(FWTrackResidualDetailView, Residuals);