MuonResiduals1DOFFitter Class Reference

#include <MuonResiduals1DOFFitter.h>

Inheritance diagram for MuonResiduals1DOFFitter:

Public Types
enum	{ kAlign = 0, kSigma, kGamma, kNPar }
enum	{ kResid = 0, kRedChi2, kNData }
Public Types inherited from MuonResidualsFitter
enum	{   kPureGaussian, kPowerLawTails, kROOTVoigt, kGaussPowerTails,   kPureGaussian2D }
enum	{   k1DOF, k5DOF, k6DOF, k6DOFrphi,   kPositionFitter, kAngleFitter, kAngleBfieldFitter }
enum	{   k1111, k1110, k1100, k1010,   k0010, k1000, k0100 }

Public Member Functions
bool	fit (Alignable *ali) override
	MuonResiduals1DOFFitter (int residualsModel, int minHits, int useResiduals, bool weightAlignment=true)
int	ndata () override
int	npar () override
double	plot (std::string name, TFileDirectory *dir, Alignable *ali) override
double	sumofweights () override
int	type () const override
Public Member Functions inherited from MuonResidualsFitter
void	computeHistogramRangeAndBinning (int which, int &nbins, double &a, double &b)
virtual void	correctBField ()=0
virtual void	correctBField (int idx_momentum, int idx_q)
TMatrixDSym	correlationMatrix ()
double	covarianceElement (int parNum1, int parNum2)
TMatrixDSym	covarianceMatrix ()
void	eraseNotSelectedResiduals ()
double	errorerror (int parNum)
void	fiducialCuts (double xMin=-80.0, double xMax=80.0, double yMin=-80.0, double yMax=80.0, bool fidcut1=false)
void	fill (double *residual)
void	fix (int parNum, bool dofix=true)
bool	fixed (int parNum)
void	histogramChi2GaussianFit (int which, double &fit_mean, double &fit_sigma)
double	loglikelihood ()
	MuonResidualsFitter (int residualsModel, int minHits, int useResiduals, bool weightAlignment=true)
int	nfixed ()
long	numResiduals () const
long	numsegments ()
int	parNum2parIdx (int parNum)
void	plotsimple (std::string name, TFileDirectory *dir, int which, double multiplier)
void	plotweighted (std::string name, TFileDirectory *dir, int which, int whichredchi2, double multiplier)
void	read (FILE *file, int which=0)
std::vector< double * >::const_iterator	residuals_begin () const
std::vector< double * >::const_iterator	residuals_end () const
int	residualsModel () const
std::vector< bool > &	selectedResidualsFlags ()
void	selectPeakResiduals (double nsigma, int nvar, int *vars)
void	selectPeakResiduals_simple (double nsigma, int nvar, int *vars)
void	setInitialValue (int parNum, double value)
void	setPrintLevel (int printLevel)
void	setStrategy (int strategy)
int	useRes (int pattern=-1)
double	value (int parNum)
void	write (FILE *file, int which=0)
virtual	~MuonResidualsFitter ()

Protected Member Functions
void	inform (TMinuit *tMinuit) override
Protected Member Functions inherited from MuonResidualsFitter
bool	dofit (void(*fcn)(int &, double *, double &, double *, int), std::vector< int > &parNum, std::vector< std::string > &parName, std::vector< double > &start, std::vector< double > &step, std::vector< double > &low, std::vector< double > &high)
void	initialize_table ()

Additional Inherited Members
Protected Attributes inherited from MuonResidualsFitter
double	m_center [20]
TMatrixDSym	m_cov
std::vector< double >	m_error
std::vector< bool >	m_fixed
double	m_loglikelihood
int	m_minHits
std::map< int, double >	m_parNum2InitValue
std::map< int, int >	m_parNum2parIdx
int	m_printLevel
double	m_radii [20]
std::vector< double * >	m_residuals
std::vector< bool >	m_residuals_ok
int	m_residualsModel
int	m_strategy
int	m_useResiduals
std::vector< double >	m_value
bool	m_weightAlignment

Detailed Description

$Date: Fri Apr 17 16:09:40 CDT 2009

Revision

1.4

Author

J. Pivarski - Texas A&M University pivar.nosp@m.ski@.nosp@m.physi.nosp@m.cs.t.nosp@m.amu.e.nosp@m.du

Definition at line 12 of file MuonResiduals1DOFFitter.h.

Member Enumeration Documentation

anonymous enum

anonymous enum

Enumerator
kAlign
kSigma
kGamma
kNPar

Definition at line 14 of file MuonResiduals1DOFFitter.h.

{ kAlign = 0, kSigma, kGamma, kNPar };

anonymous enum

anonymous enum

Enumerator
kResid
kRedChi2
kNData

Definition at line 16 of file MuonResiduals1DOFFitter.h.

{ kResid = 0, kRedChi2, kNData };

Constructor & Destructor Documentation

MuonResiduals1DOFFitter()

MuonResiduals1DOFFitter::MuonResiduals1DOFFitter ( int  residualsModel, int  minHits, int  useResiduals, bool  weightAlignment = true  )

inline

Definition at line 18 of file MuonResiduals1DOFFitter.h.

      : MuonResidualsFitter(residualsModel, minHits, useResiduals, weightAlignment) {}

Member Function Documentation

fit()

bool MuonResiduals1DOFFitter::fit ( Alignable *  ali )

overridevirtual

Implements MuonResidualsFitter.

Definition at line 64 of file MuonResiduals1DOFFitter.cc.

References MuonResidualsFitter::dofit(), MuonResidualsFitter::fixed(), LaserClient_cfi::high, MuonResidualsFitter::initialize_table(), kAlign, kGamma, MuonResidualsFitter::kGaussPowerTails, MuonResidualsFitter::kPureGaussian, kRedChi2, kResid, kSigma, LaserClient_cfi::low, MuonResidualsFitter::m_weightAlignment, MuonResiduals1DOFFitter_FCN(), N, Skims_PA_cff::name, EgammaValidation_cff::num, funct::pow(), MuonResidualsFitter::residuals_begin(), MuonResidualsFitter::residuals_end(), MuonResidualsFitter::residualsModel(), mathSSE::sqrt(), command_line::start, AlCaHLTBitMon_QueryRunRegistry::string, sumofweights(), and mps_merge::weight.

Referenced by trackingPlots.Iteration::modules().

                                                {
  initialize_table();  // if not already initialized
  sumofweights();

  double resid_sum = 0.;
  double resid_sum2 = 0.;
  double resid_N = 0.;
  int N = 0;

  for (std::vector<double *>::const_iterator resiter = residuals_begin(); resiter != residuals_end(); ++resiter) {
    const double residual = (*resiter)[MuonResiduals1DOFFitter::kResid];
    const double redchi2 = (*resiter)[MuonResiduals1DOFFitter::kRedChi2];
    double weight = 1. / redchi2;
    if (!m_weightAlignment)
      weight = 1.;

    if (!m_weightAlignment || TMath::Prob(redchi2 * 8, 8) < 0.99) {  // no spikes allowed
      if (fabs(residual) < 10.) {                                    // 10 cm
        resid_sum += weight * residual;
        resid_sum2 += weight * residual * residual;
        resid_N += weight;
        N++;
      }
    }
  }

  double resid_mean = resid_sum / resid_N;
  double resid_stdev = sqrt(resid_sum2 / resid_N - pow(resid_sum / resid_N, 2));

  resid_sum = 0.;
  resid_sum2 = 0.;
  resid_N = 0.;

  for (std::vector<double *>::const_iterator resiter = residuals_begin(); resiter != residuals_end(); ++resiter) {
    const double residual = (*resiter)[MuonResiduals1DOFFitter::kResid];
    const double redchi2 = (*resiter)[MuonResiduals1DOFFitter::kRedChi2];
    double weight = 1. / redchi2;
    if (!m_weightAlignment)
      weight = 1.;

    if (!m_weightAlignment || TMath::Prob(redchi2 * 8, 8) < 0.99) {  // no spikes allowed
      if (fabs(residual - resid_mean) < 2.5 * resid_stdev) {
        resid_sum += weight * residual;
        resid_sum2 += weight * residual * residual;
        resid_N += weight;
      }
    }
  }

  resid_mean = resid_sum / resid_N;
  resid_stdev = sqrt(resid_sum2 / resid_N - pow(resid_sum / resid_N, 2));

  std::vector<int> num;
  std::vector<std::string> name;
  std::vector<double> start;
  std::vector<double> step;
  std::vector<double> low;
  std::vector<double> high;

  if (fixed(kAlign)) {
    num.push_back(kAlign);
    name.push_back(std::string("Align"));
    start.push_back(0.);
    step.push_back(0.01 * resid_stdev);
    low.push_back(0.);
    high.push_back(0.);
  } else {
    num.push_back(kAlign);
    name.push_back(std::string("Align"));
    start.push_back(resid_mean);
    step.push_back(0.01 * resid_stdev);
    low.push_back(0.);
    high.push_back(0.);
  }
  num.push_back(kSigma);
  name.push_back(std::string("Sigma"));
  start.push_back(resid_stdev);
  step.push_back(0.01 * resid_stdev);
  low.push_back(0.);
  high.push_back(0.);
  if (residualsModel() != kPureGaussian && residualsModel() != kGaussPowerTails) {
    num.push_back(kGamma);
    name.push_back(std::string("Gamma"));
    start.push_back(0.1 * resid_stdev);
    step.push_back(0.01 * resid_stdev);
    low.push_back(0.);
    high.push_back(0.);
  }

  return dofit(&MuonResiduals1DOFFitter_FCN, num, name, start, step, low, high);
}

inform()

void MuonResiduals1DOFFitter::inform ( TMinuit *  tMinuit )

overrideprotectedvirtual

Implements MuonResidualsFitter.

Definition at line 9 of file MuonResiduals1DOFFitter.cc.

References MuonResiduals1DOFFitter_TMinuit.

{ MuonResiduals1DOFFitter_TMinuit = tMinuit; }

ndata()

int MuonResiduals1DOFFitter::ndata ( )

inlineoverridevirtual

Implements MuonResidualsFitter.

Definition at line 33 of file MuonResiduals1DOFFitter.h.

References kNData.

{ return kNData; }

npar()

int MuonResiduals1DOFFitter::npar ( )

inlineoverridevirtual

Implements MuonResidualsFitter.

Definition at line 23 of file MuonResiduals1DOFFitter.h.

References cms::cuda::assert(), MuonResidualsFitter::kGaussPowerTails, kNPar, MuonResidualsFitter::kPowerLawTails, MuonResidualsFitter::kPureGaussian, MuonResidualsFitter::kROOTVoigt, and MuonResidualsFitter::residualsModel().

Referenced by plot().

                      {
    if (residualsModel() == kPureGaussian || residualsModel() == kGaussPowerTails)
      return kNPar - 1;
    else if (residualsModel() == kPowerLawTails)
      return kNPar;
    else if (residualsModel() == kROOTVoigt)
      return kNPar;
    else
      assert(false);
  }

plot()

double MuonResiduals1DOFFitter::plot ( std::string  name, TFileDirectory *  dir, Alignable *  ali  )

overridevirtual

Implements MuonResidualsFitter.

Definition at line 156 of file MuonResiduals1DOFFitter.cc.

References cms::cuda::assert(), hltPixelTracks_cff::chi2, DeadROC_duringRun::dir, mps_fire::i, kAlign, kGamma, MuonResidualsFitter::kGaussPowerTails, MuonResidualsFitter::kPowerLawTails, MuonResidualsFitter::kPureGaussian, kRedChi2, kResid, MuonResidualsFitter::kROOTVoigt, kSigma, MuonResidualsFitter::m_weightAlignment, MuonResiduals1DOFFitter_sum_of_weights, MuonResidualsFitter_GaussPowerTails_TF1(), MuonResidualsFitter_powerLawTails_TF1(), MuonResidualsFitter_pureGaussian_TF1(), MuonResidualsFitter_ROOTVoigt_TF1(), Skims_PA_cff::name, ndof, npar(), funct::pow(), MuonResidualsFitter::residuals_begin(), MuonResidualsFitter::residuals_end(), MuonResidualsFitter::residualsModel(), sumofweights(), MuonResidualsFitter::value(), and protons_cff::xi.

                                                                                        {
  sumofweights();

  std::stringstream name_residual, name_residual_raw;
  name_residual << name << "_residual";
  name_residual_raw << name << "_residual_raw";

  double min_residual = -100.;
  double max_residual = 100.;
  TH1F *hist_residual = dir->make<TH1F>(name_residual.str().c_str(), "", 100, min_residual, max_residual);
  TH1F *hist_residual_raw = dir->make<TH1F>(name_residual_raw.str().c_str(), "", 100, min_residual, max_residual);

  name_residual << "_fit";
  TF1 *fit_residual = nullptr;
  if (residualsModel() == kPureGaussian) {
    fit_residual =
        new TF1(name_residual.str().c_str(), MuonResidualsFitter_pureGaussian_TF1, min_residual, max_residual, 3);
    fit_residual->SetParameters(MuonResiduals1DOFFitter_sum_of_weights * (max_residual - min_residual) / 100.,
                                10. * value(kAlign),
                                10. * value(kSigma));
  } else if (residualsModel() == kPowerLawTails) {
    fit_residual =
        new TF1(name_residual.str().c_str(), MuonResidualsFitter_powerLawTails_TF1, min_residual, max_residual, 4);
    fit_residual->SetParameters(MuonResiduals1DOFFitter_sum_of_weights * (max_residual - min_residual) / 100.,
                                10. * value(kAlign),
                                10. * value(kSigma),
                                10. * value(kGamma));
  } else if (residualsModel() == kROOTVoigt) {
    fit_residual =
        new TF1(name_residual.str().c_str(), MuonResidualsFitter_ROOTVoigt_TF1, min_residual, max_residual, 4);
    fit_residual->SetParameters(MuonResiduals1DOFFitter_sum_of_weights * (max_residual - min_residual) / 100.,
                                10. * value(kAlign),
                                10. * value(kSigma),
                                10. * value(kGamma));
  } else if (residualsModel() == kGaussPowerTails) {
    fit_residual =
        new TF1(name_residual.str().c_str(), MuonResidualsFitter_GaussPowerTails_TF1, min_residual, max_residual, 3);
    fit_residual->SetParameters(MuonResiduals1DOFFitter_sum_of_weights * (max_residual - min_residual) / 100.,
                                10. * value(kAlign),
                                10. * value(kSigma));
  } else {
    assert(false);
  }

  fit_residual->SetLineColor(2);
  fit_residual->SetLineWidth(2);
  fit_residual->Write();

  for (std::vector<double *>::const_iterator resiter = residuals_begin(); resiter != residuals_end(); ++resiter) {
    const double resid = (*resiter)[MuonResiduals1DOFFitter::kResid];
    const double redchi2 = (*resiter)[MuonResiduals1DOFFitter::kRedChi2];
    double weight = 1. / redchi2;
    if (!m_weightAlignment)
      weight = 1.;

    if (!m_weightAlignment || TMath::Prob(redchi2 * 8, 8) < 0.99) {  // no spikes allowed
      hist_residual->Fill(10. * (resid + value(kAlign)), weight);
    }

    hist_residual_raw->Fill(10. * resid);
  }

  double chi2 = 0.;
  double ndof = 0.;
  for (int i = 1; i <= hist_residual->GetNbinsX(); i++) {
    double xi = hist_residual->GetBinCenter(i);
    double yi = hist_residual->GetBinContent(i);
    double yerri = hist_residual->GetBinError(i);
    double yth = fit_residual->Eval(xi);
    if (yerri > 0.) {
      chi2 += pow((yth - yi) / yerri, 2);
      ndof += 1.;
    }
  }
  ndof -= npar();

  return (ndof > 0. ? chi2 / ndof : -1.);
}

sumofweights()

double MuonResiduals1DOFFitter::sumofweights ( )

overridevirtual

Implements MuonResidualsFitter.

Definition at line 45 of file MuonResiduals1DOFFitter.cc.

References kRedChi2, MuonResidualsFitter::m_weightAlignment, MuonResiduals1DOFFitter_number_of_hits, MuonResiduals1DOFFitter_sum_of_weights, MuonResiduals1DOFFitter_weightAlignment, MuonResidualsFitter::residuals_begin(), and MuonResidualsFitter::residuals_end().

Referenced by fit(), and plot().

                                             {
  MuonResiduals1DOFFitter_sum_of_weights = 0.;
  MuonResiduals1DOFFitter_number_of_hits = 0.;
  MuonResiduals1DOFFitter_weightAlignment = m_weightAlignment;
  for (std::vector<double *>::const_iterator resiter = residuals_begin(); resiter != residuals_end(); ++resiter) {
    if (m_weightAlignment) {
      double redchi2 = (*resiter)[MuonResiduals1DOFFitter::kRedChi2];
      if (TMath::Prob(redchi2 * 8, 8) < 0.99) {  // no spikes allowed
        MuonResiduals1DOFFitter_sum_of_weights += 1. / redchi2;
        MuonResiduals1DOFFitter_number_of_hits += 1.;
      }
    } else {
      MuonResiduals1DOFFitter_sum_of_weights += 1.;
      MuonResiduals1DOFFitter_number_of_hits += 1.;
    }
  }
  return MuonResiduals1DOFFitter_sum_of_weights;
}

type()

int MuonResiduals1DOFFitter::type ( ) const

inlineoverridevirtual

Implements MuonResidualsFitter.

Definition at line 21 of file MuonResiduals1DOFFitter.h.

References MuonResidualsFitter::k1DOF.

{ return MuonResidualsFitter::k1DOF; }