CMS 3D CMS Logo

Functions

/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.cc File Reference

#include "DQM/BeamMonitor/plugins/Vx3DHLTAnalyzer.h"
#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/Framework/interface/LuminosityBlock.h"
#include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHitCollection.h"
#include "DataFormats/TrackReco/interface/Track.h"
#include "DataFormats/TrackReco/interface/TrackFwd.h"
#include "DataFormats/VertexReco/interface/VertexFwd.h"
#include "DataFormats/VertexReco/interface/Vertex.h"
#include <TFitterMinuit.h>

Go to the source code of this file.

Functions

 DEFINE_FWK_MODULE (Vx3DHLTAnalyzer)
void Gauss3DFunc (int &, double *, double &fval, double *par, int)

Function Documentation

DEFINE_FWK_MODULE ( Vx3DHLTAnalyzer  )
void Gauss3DFunc ( int &  ,
double *  ,
double &  fval,
double *  par,
int   
)

Definition at line 176 of file Vx3DHLTAnalyzer.cc.

References considerVxCovariance, counterVx, DIM, i, funct::log(), maxLongLength, maxTransRadius, pi, mathSSE::sqrt(), Vertices, VxErrCorr, ExpressReco_HICollisions_FallBack::x, xPos, ExpressReco_HICollisions_FallBack::y, yPos, z, and zPos.

Referenced by Vx3DHLTAnalyzer::MyFit().

{
  double K[DIM][DIM]; // Covariance Matrix
  double M[DIM][DIM]; // K^-1
  double det;
  double sumlog = 0.;

//   par[0] = K(0,0) --> Var[X]
//   par[1] = K(1,1) --> Var[Y]
//   par[2] = K(2,2) --> Var[Z]
//   par[3] = K(0,1) = K(1,0) --> Cov[X,Y]
//   par[4] = K(1,2) = K(2,1) --> Cov[Y,Z] --> dy/dz
//   par[5] = K(0,2) = K(2,0) --> Cov[X,Z] --> dx/dz
//   par[6] = mean x
//   par[7] = mean y
//   par[8] = mean z

  counterVx = 0;
  for (unsigned int i = 0; i < Vertices.size(); i++)
    {
      if ((std::sqrt((Vertices[i].x-xPos)*(Vertices[i].x-xPos) + (Vertices[i].y-yPos)*(Vertices[i].y-yPos)) <= maxTransRadius) &&
          (std::fabs(Vertices[i].z-zPos) <= maxLongLength))
        {
          if (considerVxCovariance == true)
            {
              K[0][0] = std::fabs(par[0]) + VxErrCorr*VxErrCorr * std::fabs(Vertices[i].Covariance[0][0]);
              K[1][1] = std::fabs(par[1]) + VxErrCorr*VxErrCorr * std::fabs(Vertices[i].Covariance[1][1]);
              K[2][2] = std::fabs(par[2]) + VxErrCorr*VxErrCorr * std::fabs(Vertices[i].Covariance[2][2]);
              K[0][1] = K[1][0] = par[3] + VxErrCorr*VxErrCorr * Vertices[i].Covariance[0][1];
              K[1][2] = K[2][1] = par[4]*(std::fabs(par[2])-std::fabs(par[1])) - par[5]*par[3] + VxErrCorr*VxErrCorr * Vertices[i].Covariance[1][2];
              K[0][2] = K[2][0] = par[5]*(std::fabs(par[2])-std::fabs(par[0])) - par[4]*par[3] + VxErrCorr*VxErrCorr * Vertices[i].Covariance[0][2];
            }
          else
            {
              K[0][0] = std::fabs(par[0]);
              K[1][1] = std::fabs(par[1]);
              K[2][2] = std::fabs(par[2]);
              K[0][1] = K[1][0] = par[3];
              K[1][2] = K[2][1] = par[4]*(std::fabs(par[2])-std::fabs(par[1])) - par[5]*par[3];
              K[0][2] = K[2][0] = par[5]*(std::fabs(par[2])-std::fabs(par[0])) - par[4]*par[3];
            }

          det = K[0][0]*(K[1][1]*K[2][2] - K[1][2]*K[1][2]) -
                K[0][1]*(K[0][1]*K[2][2] - K[0][2]*K[1][2]) +
                K[0][2]*(K[0][1]*K[1][2] - K[0][2]*K[1][1]);

          M[0][0] = (K[1][1]*K[2][2] - K[1][2]*K[1][2]) / det;
          M[1][1] = (K[0][0]*K[2][2] - K[0][2]*K[0][2]) / det;
          M[2][2] = (K[0][0]*K[1][1] - K[0][1]*K[0][1]) / det;
          M[0][1] = M[1][0] = (K[0][2]*K[1][2] - K[0][1]*K[2][2]) / det;
          M[1][2] = M[2][1] = (K[0][2]*K[0][1] - K[1][2]*K[0][0]) / det;
          M[0][2] = M[2][0] = (K[0][1]*K[1][2] - K[0][2]*K[1][1]) / det;
          
          sumlog += double(DIM)*std::log(2.*pi) + std::log(std::fabs(det)) +
            (M[0][0]*(Vertices[i].x-par[6])*(Vertices[i].x-par[6]) +
             M[1][1]*(Vertices[i].y-par[7])*(Vertices[i].y-par[7]) +
             M[2][2]*(Vertices[i].z-par[8])*(Vertices[i].z-par[8]) +
             2.*M[0][1]*(Vertices[i].x-par[6])*(Vertices[i].y-par[7]) +
             2.*M[1][2]*(Vertices[i].y-par[7])*(Vertices[i].z-par[8]) +
             2.*M[0][2]*(Vertices[i].x-par[6])*(Vertices[i].z-par[8]));
          
          counterVx++;
        }
    }
  
  fval = sumlog;
}