CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_6/src/DataFormats/TrackReco/src/TrackResiduals.cc

Go to the documentation of this file.
00001 #include <math.h>
00002 #include "DataFormats/TrackReco/interface/HitPattern.h"
00003 #include "DataFormats/TrackReco/interface/TrackResiduals.h"
00004 #include <cstdio>
00005 #include <cstring>
00006 
00007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00008 
00009 using namespace reco;
00010 
00011 TrackResiduals::TrackResiduals () : residualType(X_Y_RESIDUALS)
00012 {
00013      memset(residuals_, 0, sizeof(residuals_));
00014 }
00015 
00016 TrackResiduals::TrackResiduals (enum ResidualType type) : residualType(type)
00017 {
00018      memset(residuals_, 0, sizeof(residuals_));
00019 }
00020 
00021 void TrackResiduals::setResidualType (enum ResidualType type)
00022 {
00023      residualType = type;
00024 }
00025 
00026 void TrackResiduals::setResidualXY (int idx, double residualX, double residualY)
00027 {
00028      assert(residualType == X_Y_RESIDUALS);
00029      if (idx>=numResiduals) {
00030        edm::LogWarning("TrackResiduals")<<" setting residual over the array size.";
00031        return;}
00032      residuals_[idx] = (pack_residual(residualX) << 4) | pack_residual(residualY);
00033 }
00034 
00035 double TrackResiduals::residualX (int i) const
00036 {
00037      switch (residualType) {
00038      case X_Y_RESIDUALS:
00039           return unpack_residual(residuals_[i] >> 4);
00040      case X_Y_PULLS:
00041           return unpack_pull(residuals_[i] >> 4);
00042      default:
00043           assert(0);
00044      }
00045      return 0;
00046 }
00047 
00048 double TrackResiduals::residualY (int i) const
00049 {
00050      switch (residualType) {
00051      case X_Y_RESIDUALS:
00052           return unpack_residual(residuals_[i] & 0x0f);
00053      case X_Y_PULLS:
00054           return unpack_pull(residuals_[i] & 0x0f);
00055      default:
00056           assert(0);
00057      }
00058      return 0;
00059 }
00060 
00061 static int index_to_hitpattern (int i_hitpattern, const HitPattern &h)
00062 {
00063      int i_residuals = 0;
00064      assert(i_hitpattern < h.numberOfHits());
00065      if (!h.validHitFilter(h.getHitPattern(i_hitpattern))) 
00066           // asking for residual of invalid hit...
00067           return -999;
00068      for (int i = 0; i <= i_hitpattern; i++) {
00069           if (h.validHitFilter(h.getHitPattern(i)))
00070                i_residuals++;
00071      }
00072      assert(i_residuals > 0);
00073      return i_residuals - 1;
00074 }
00075 
00076 double TrackResiduals::residualX (int i, const HitPattern &h) const
00077 {
00078      int idx = index_to_hitpattern(i, h);
00079      if (idx == -999)
00080           return -999;
00081      return residualX(idx);
00082 }
00083 
00084 double TrackResiduals::residualY (int i, const HitPattern &h) const
00085 {
00086      int idx = index_to_hitpattern(i, h);
00087      if (idx == -999)
00088           return -999;
00089      return residualY(idx);
00090 }
00091 
00092 void TrackResiduals::setPullXY (int idx, double pullX, double pullY)
00093 {
00094      assert(residualType == X_Y_PULLS);
00095      if (idx>=numResiduals) {
00096        edm::LogWarning("TrackResiduals")<<" setting pulls over the array size.";
00097        return;}
00098 
00099      residuals_[idx] = (pack_pull(pullX) << 4) | pack_pull(pullY);
00100 }
00101 
00102 static const double pull_char_to_double[8][2] = { 
00103      { 0,   0.5 },
00104      { 0.5, 1   },
00105      { 1,   1.5 },
00106      { 1.5, 2   },
00107      { 2,   2.5 },
00108      { 2.5, 3.5 },
00109      { 3.5, 4.5 },
00110      { 4.5, 5.5 },
00111 };
00112 
00113 double TrackResiduals::unpack_pull (unsigned char pull)
00114 {
00115      int sgn = 1 - 2 * ((pull & 0x08) >> 3);
00116      unsigned char mag = pull & 0x07;
00117      return sgn * 
00118           (pull_char_to_double[mag][0] + pull_char_to_double[mag][1]) / 2;
00119 }
00120 
00121 unsigned char TrackResiduals::pack_pull (double pull)
00122 {
00123      unsigned char sgn = (pull < 0) * 0x08; // 1xxx is -abs(0xxx)
00124      int mag = -1;
00125      while (++mag < 8 && pull_char_to_double[mag][1] < fabs(pull));
00126      return sgn + mag;
00127 }
00128 
00129 static const double residual_char_to_double[8][2] = { 
00130      { 0,   0.5 },
00131      { 0.5, 1   },
00132      { 1,   1.5 },
00133      { 1.5, 2   },
00134      { 2,   2.5 },
00135      { 2.5, 3.5 },
00136      { 3.5, 4.5 },
00137      { 4.5, 5.5 },
00138 };
00139 
00140 double TrackResiduals::unpack_residual (unsigned char pull)
00141 {
00142      int sgn = 1 - 2 * ((pull & 0x08) >> 3);
00143      unsigned char mag = pull & 0x07;
00144      return sgn * 
00145           (pull_char_to_double[mag][0] + pull_char_to_double[mag][1]) / 2;
00146 }
00147 
00148 unsigned char TrackResiduals::pack_residual (double pull)
00149 {
00150      unsigned char sgn = (pull < 0) * 0x08; // 1xxx is -abs(0xxx)
00151      int mag = -1;
00152      while (++mag < 8 && pull_char_to_double[mag][1] < fabs(pull));
00153      return sgn + mag;
00154 }
00155 
00156 void TrackResiduals::print (std::ostream &stream) const
00157 {
00158      stream << "TrackResiduals" << std::endl;
00159      std::ios_base::fmtflags flags = stream.flags();
00160      stream.setf ( std::ios_base::hex, std::ios_base::basefield );  
00161      stream.setf ( std::ios_base::showbase );               
00162      for (int i = 0; i < numResiduals; i++) {
00163           unsigned char residual = residuals_[i];
00164           printf("0x%x\n", residual);
00165 //        stream << residual << std::endl;
00166      }
00167      stream.flags(flags);
00168 }
00169 
00170 void TrackResiduals::print (const HitPattern &h, std::ostream &stream) const
00171 {
00172      stream << "TrackResiduals" << std::endl;
00173      for (int i = 0; i < h.numberOfHits(); i++) {
00174           stream << (h.validHitFilter(h.getHitPattern(i)) ? 
00175                      "valid hit:   " : "invalid hit: ") << 
00176                "( " << residualX(i, h) << " , " << residualY(i, h) << " )" << 
00177                std::endl;
00178      }
00179 }