QratCscClusterFitter.cxx File Reference

#include "QratCscClusterFitter.h"
#include <iomanip>
#include <sstream>
#include "EventPrimitives/EventPrimitives.h"
#include "EventPrimitives/EventPrimitivesHelpers.h"
#include "MuonPrepRawData/CscClusterStatus.h"
#include "MuonPrepRawData/CscPrepData.h"
#include "MuonPrepRawData/CscStripPrepData.h"
#include "MuonReadoutGeometry/CscReadoutElement.h"
#include "TrkEventPrimitives/LocalDirection.h"
#include "TrkEventPrimitives/ParamDefs.h"
#include "CxxUtils/trapping_fp.h"

Go to the source code of this file.

Typedefs
using	Results = std::vector<Result>

Enumerations
enum	CscStation { UNKNOWN_STATION , CSS , CSL }
enum	CscPlane {   CSS_ETA , CSL_ETA , CSS_PHI , CSL_PHI ,   UNKNOWN_PLANE }
enum	CscClusterStatus
	Enum to represent the cluster status - see the specific enum values for more details. More...

Functions
int	qrat_correction (CscPlane plane, double qrat, double &cor, double &dcordqrat)
int	qrat_interpolation (double qrmin, const std::vector< double > &corvals, double qrat, double &cor, double &dcordqrat)
int	qrat_atanh (const double a, const double b, double c, const double x0, double qrat, double &cor, double &dcordqrat)
	Calculate QRAT correction from inverse hyperbolic tangent based on a fit to the plot of pos x vs charge ratio qrat to the function qrat(x) = a + b* tanh(c*(x-x0)) Here we use the inverse of the above function x = atanh((qrat-a)/b)/c+x0 to obtain positions from charge ratios.

Typedef Documentation

Results

using Results = std::vector<Result>

Definition at line 27 of file QratCscClusterFitter.cxx.

Enumeration Type Documentation

CscClusterStatus

enum Muon::CscClusterStatus

Enum to represent the cluster status - see the specific enum values for more details.

Definition at line 23 of file CscClusterStatus.h.

                        {

    CscStatusUnspoiled =0,

    CscStatusSimple =1,
 
    // Spoiled cluster list, it should be handled with care in reconstruction
    
    CscStatusEdge =2,

    CscStatusMultiPeak =3,

    CscStatusNarrow =4,

    CscStatusWide =5,

    CscStatusSkewed =6,

    CscStatusQratInconsistent =7,
 
    // Parabola interpolation failed in peak strip
    // OR either left or right strip is hot or dead
    CscStatusStripFitFailed = 8,
 
    // any of left/peak/right strip is saturated...
    CscStatusSaturated = 9,

    CscStatusSplitUnspoiled =10,

    CscStatusSplitSimple =11,
 
    // Spoiled cluster list, it should be handled with care in reconstruction
    
    CscStatusSplitEdge =12,
    
    CscStatusSplitMultiPeak =13,

    CscStatusSplitNarrow =14,

    CscStatusSplitWide =15,

    CscStatusSplitSkewed =16,
    
    CscStatusSplitQratInconsistent =17,
    
    // Parabola interpolation failed in peak strip
    // OR either left or right strip is hot or dead
    CscStatusSplitStripFitFailed = 18,
 
    // any of left/peak/right strip is saturated...
    CscStatusSplitSaturated = 19,

    CscStatusUndefined =255
  };

CscPlane

enum CscPlane

Enumerator
CSS_ETA
CSL_ETA
CSS_PHI
CSL_PHI
UNKNOWN_PLANE

Definition at line 30 of file QratCscClusterFitter.cxx.

{ CSS_ETA, CSL_ETA, CSS_PHI, CSL_PHI, UNKNOWN_PLANE };

CscStation

enum CscStation

Enumerator
UNKNOWN_STATION
CSS
CSL

Definition at line 29 of file QratCscClusterFitter.cxx.

{ UNKNOWN_STATION, CSS, CSL };

Function Documentation

qrat_atanh()

int qrat_atanh	(	const double	a,
		const double	b,
		double	c,
		const double	x0,
		double	qrat,
		double &	cor,
		double &	dcordqrat )

Calculate QRAT correction from inverse hyperbolic tangent based on a fit to the plot of pos x vs charge ratio qrat to the function qrat(x) = a + b* tanh(c*(x-x0)) Here we use the inverse of the above function x = atanh((qrat-a)/b)/c+x0 to obtain positions from charge ratios.

Output is shifted by 0.5, so the center of the strip is at x=0.

Parameters

a	= parameter of correction function
b	= parameter of correction function
c	= parameter of correction function
x0	= parameter of correction function
qrat	= input value of qrat = Q_right/Q_peak
cor	= output of corrected position in the range -0.5 to 0.5
dcordqrat	= derivative of cor w.r.t. qrat for error estimates

Returns

0 if no error, 1 if qrat is negative

Definition at line 194 of file QratCscClusterFitter.cxx.

                                                                                                                       {
    if (qrat <= 0) return 1;  // avoid trouble in error calculation
    // minimum qrat value (at pos=-1) or use pos = -0.5?
    double qrmin = a + b * tanh(c * (-1 - x0));
 
    // Treat any QRAT below the minimum as if it were at the minumum.
    if (qrat < qrmin) qrat = qrmin;
    double z = (qrat - a) / b;
    cor = atanh(z) / c + x0;
    dcordqrat = 1.0 / (1.0 - z * z) / b / c;
 
    return 0;
}

qrat_correction()

int qrat_correction	(	CscPlane	plane,
		double	qrat,
		double &	cor,
		double &	dcordqrat )

Definition at line 62 of file QratCscClusterFitter.cxx.

                                                                                 {
    std::vector<double> pfac;
    if (plane == CSS_ETA) {
        if (qrat < 0.095) return 1;
        if (qrat > 1.01) return 1;
        if (qrat < 0.15) {
            pfac.push_back(-5.483);
            pfac.push_back(130.6);
            pfac.push_back(-1296.);
            pfac.push_back(5994.);
            pfac.push_back(-10580.);
        } else if (qrat < 0.30) {
            pfac.push_back(-0.9225);
            pfac.push_back(5.569);
            pfac.push_back(2.908);
            pfac.push_back(-66.19);
            pfac.push_back(108.5);
        } else {
            pfac.push_back(-0.4246);
            pfac.push_back(2.619);
            pfac.push_back(-3.642);
            pfac.push_back(2.837);
            pfac.push_back(-0.8916);
        }
    } else if (plane == CSL_ETA) {
        if (qrat < 0.105) return 1;
        if (qrat > 1.01) return 1;
        if (qrat < 0.25) {
            pfac.push_back(-2.823);
            pfac.push_back(42.71);
            pfac.push_back(-279.1);
            pfac.push_back(879.2);
            pfac.push_back(-1062.);
        } else {
            pfac.push_back(-0.5409);
            pfac.push_back(3.110);
            pfac.push_back(-4.630);
            pfac.push_back(3.800);
            pfac.push_back(-1.241);
        }
    } else {
        return 3;
    }
    cor = 0.0;
    dcordqrat = 0.0;
    double term = 1.0;
    for (unsigned int ipow = 0; ipow < pfac.size(); ++ipow) {
        dcordqrat += ipow * pfac[ipow] * term;
        if (ipow) term *= qrat;
        cor += pfac[ipow] * term;
    }
    return 0;
}

qrat_interpolation()

int qrat_interpolation	(	double	qrmin,
		const std::vector< double > &	corvals,
		double	qrat,
		double &	cor,
		double &	dcordqrat )

Definition at line 128 of file QratCscClusterFitter.cxx.

                                                                                                                    {
    int nbin = corvals.size();
    if (!nbin) return 1;
    // Treat any QRAT below the minimum as if it were at the minumum.
    if (qrat < qrmin) qrat = qrmin;
    // Find the bin holding dqrat.
    double dqrat = 1.0 / nbin;
    int bin = int(qrat / dqrat);
    // Extract the value for this bin (x1) and the preceding (x0)
    // and the following (x2).
    double x1 = corvals[bin];
    if (x1 == 0.0) return 2;
    double x0 = bin == 0 ? 0.0 : corvals[bin - 1];
    double x2 = bin == nbin - 1 ? 1.0 : corvals[bin + 1];
    // Assign the qrat values for these bins. If this is the first
    // bin with data, then use the input qrmin.
    double qrat0 = dqrat * bin;
    double qrat1 = qrat0 + dqrat;
 
    if (x0 == 0.0) {
        if (qrmin > qrat0 && qrmin < qrat1) { qrat0 = qrmin; }
    }
    // Calculate correction and derivative.
    // Use quadratic interpolation with the high edge of this bin,
    // the preceding and the following. For the last bin, use linear
    // interpolation.
    if (bin == nbin - 1) {
        double a = x0;
        double b = (x1 - x0) / dqrat;
        cor = a + b * (qrat - qrat0);
        dcordqrat = b;
    } else {
        double d0 = qrat1 - qrat0;
        double d = dqrat;
        double w = 1.0 / (d0 * d * d + d0 * d0 * d);
        double a = x1;
        double b = w * d * d * (x1 - x0) + w * d0 * d0 * (x2 - x1);
        double c = w * d0 * (x2 - x1) - w * d * (x1 - x0);
        cor = a + b * (qrat - qrat1) + c * (qrat - qrat1) * (qrat - qrat1);
        dcordqrat = b + 2.0 * c * (qrat - qrat1);
    }
    if (cor < 0.0) cor = 0.0;
    // Shift correction to center of bin.
    cor -= 0.50;
    return 0;
}