#include <string.h>
#include "AliMultiplicity.h"
#include "AliLog.h"
+#include "AliRefArray.h"
ClassImp(AliMultiplicity)
TObject(),
fNtracks(0),
fNsingle(0),
+//
+ fDPhiWindow2(0.08*0.08),
+ fDThetaWindow2(0.025*0.025),
+ fDPhiShift(0.0045),
+ fNStdDev(1.0),
+//
fLabels(0),
fLabelsL2(0),
- fUsedClusS(0),
+ fUsedClusS(0),
fUsedClusT(0),
fTh(0),
fPhi(0),
// Default Constructor
fFiredChips[0] = 0;
fFiredChips[1] = 0;
+ for (int il=2;il--;) fSCl2Tracks[il] = fTCl2Tracks[il][0] = fTCl2Tracks[il][1] = 0;
for(Int_t ilayer = 0; ilayer < 6; ilayer++)fITSClusters[ilayer] = 0;
}
TObject(),
fNtracks(ntr),
fNsingle(ns),
+ //
+ fDPhiWindow2(0.08*0.08),
+ fDThetaWindow2(0.025*0.025),
+ fDPhiShift(0.0045),
+ fNStdDev(1.0),
+ //
fLabels(0),
fLabelsL2(0),
fUsedClusS(0),
fClusterFiredChips(1200)
{
// Standard constructor
+ for (int il=2;il--;) fSCl2Tracks[il] = fTCl2Tracks[il][0] = fTCl2Tracks[il][1] = 0;
if(ntr>0){
fLabels = new Int_t[ntr];
fLabelsL2 = new Int_t[ntr];
TObject(),
fNtracks(ntr),
fNsingle(ns),
+ //
+ fDPhiWindow2(0.08*0.08),
+ fDThetaWindow2(0.025*0.025),
+ fDPhiShift(0.0045),
+ fNStdDev(1.0),
+ //
fLabels(0),
fLabelsL2(0),
fUsedClusS(0),
fClusterFiredChips(1200)
{
// Standard constructor to create the arrays w/o filling
+ for (int il=2;il--;) fSCl2Tracks[il] = fTCl2Tracks[il][0] = fTCl2Tracks[il][1] = 0;
if(ntr>0){
fLabels = new Int_t[ntr];
fLabelsL2 = new Int_t[ntr];
fPhi = new Double_t [ntr];
fDeltTh = new Double_t [ntr];
fDeltPhi = new Double_t [ntr];
- fUsedClusT = new ULong64_t[ntr];
-
for(Int_t i=fNtracks;i--;){
fTh[i]=fPhi[i]=fDeltTh[i]=fDeltPhi[i] = 0;
fLabels[i] = fLabelsL2[i] = 0;
fThsingle = new Double_t [ns];
fPhisingle = new Double_t [ns];
fLabelssingle = new Int_t [ns];
- fUsedClusS = new UInt_t[ns];
for(Int_t i=fNsingle;i--;) fThsingle[i] = fPhisingle[i] = fLabelssingle[i] = 0;
}
fFiredChips[0] = nfcL1;
TObject(m),
fNtracks(m.fNtracks),
fNsingle(m.fNsingle),
+ //
+ fDPhiWindow2(0.08*0.08),
+ fDThetaWindow2(0.025*0.025),
+ fDPhiShift(0.0045),
+ fNStdDev(1.0),
+ //
fLabels(0),
fLabelsL2(0),
fUsedClusS(0),
fClusterFiredChips(1200)
{
// copy constructor
+ for (int il=2;il--;) fSCl2Tracks[il] = fTCl2Tracks[il][0] = fTCl2Tracks[il][1] = 0;
Duplicate(m);
}
if(fLabelssingle)delete [] fLabelssingle;fLabelssingle = 0;
if(fUsedClusS) delete[] fUsedClusS; fUsedClusS = 0;
if(fUsedClusT) delete[] fUsedClusT; fUsedClusT = 0;
-
+ for (int il=2;il--;) {
+ if (fSCl2Tracks[il]) delete fSCl2Tracks[il]; fSCl2Tracks[il] = 0;
+ if (fTCl2Tracks[il][0]) delete fTCl2Tracks[il][0]; fTCl2Tracks[il][0] = 0;
+ if (fTCl2Tracks[il][1]) delete fTCl2Tracks[il][1]; fTCl2Tracks[il][1] = 0;
+ }
Duplicate(m);
-
+ //
return *this;
}
fDeltPhi = new Double_t[fNtracks];
fLabels = new Int_t[fNtracks];
fLabelsL2 = new Int_t[fNtracks];
- if (m.fUsedClusT) fUsedClusT = new ULong64_t[fNtracks];
- else fUsedClusT = 0;
+ if (m.fUsedClusT) fUsedClusT = new ULong64_t[fNtracks]; else fUsedClusT = 0;
if(m.fTh)memcpy(fTh,m.fTh,fNtracks*sizeof(Double_t));
if(m.fPhi)memcpy(fPhi,m.fPhi,fNtracks*sizeof(Double_t));
if(m.fDeltTh)memcpy(fDeltTh,m.fDeltTh,fNtracks*sizeof(Double_t));
if(m.fLabels)memcpy(fLabels,m.fLabels,fNtracks*sizeof(Int_t));
if(m.fLabelsL2)memcpy(fLabelsL2,m.fLabelsL2,fNtracks*sizeof(Int_t));
if(fUsedClusT) memcpy(fUsedClusT,m.fUsedClusT,fNtracks*sizeof(ULong64_t));
+ for (int i=2;i--;) for (int j=2;j--;) if (m.fTCl2Tracks[i][j]) fTCl2Tracks[i][j] = new AliRefArray(*m.fTCl2Tracks[i][j]);
}
else {
fTh = 0;
if(m.fPhisingle)memcpy(fPhisingle,m.fPhisingle,fNsingle*sizeof(Double_t));
if(m.fLabelssingle)memcpy(fLabelssingle,m.fLabelssingle,fNsingle*sizeof(Int_t));
if(fUsedClusS) memcpy(fUsedClusS,m.fUsedClusS,fNsingle*sizeof(UInt_t));
+ for (int i=2;i--;) if (m.fSCl2Tracks[i]) fSCl2Tracks[i] = new AliRefArray(*m.fSCl2Tracks[i]);
}
else {
fThsingle = 0;
for(Int_t ilayer = 0; ilayer < 6; ilayer++){
fITSClusters[ilayer] = m.fITSClusters[ilayer];
}
-
+ fDPhiWindow2 = m.fDPhiWindow2;
+ fDThetaWindow2 = m.fDThetaWindow2;
+ fDPhiShift = fDPhiShift;
+ fNStdDev = fNStdDev;
fFastOrFiredChips = m.fFastOrFiredChips;
fClusterFiredChips = m.fClusterFiredChips;
}
if(fLabelssingle)delete [] fLabelssingle;fLabelssingle = 0;
if(fUsedClusS) delete[] fUsedClusS; fUsedClusS = 0;
if(fUsedClusT) delete[] fUsedClusT; fUsedClusT = 0;
+ for (int il=2;il--;) {
+ if (fSCl2Tracks[il]) delete fSCl2Tracks[il]; fSCl2Tracks[il] = 0;
+ if (fTCl2Tracks[il][0]) delete fTCl2Tracks[il][0]; fTCl2Tracks[il][0] = 0;
+ if (fTCl2Tracks[il][1]) delete fTCl2Tracks[il][1]; fTCl2Tracks[il][1] = 0;
+ }
}
//______________________________________________________________________
if(fLabelssingle)delete [] fLabelssingle;fLabelssingle = 0;
if(fUsedClusS) delete[] fUsedClusS; fUsedClusS = 0;
if(fUsedClusT) delete[] fUsedClusT; fUsedClusT = 0;
+ for (int il=2;il--;) {
+ if (fSCl2Tracks[il]) delete fSCl2Tracks[il]; fSCl2Tracks[il] = 0;
+ if (fTCl2Tracks[il][0]) delete fTCl2Tracks[il][0]; fTCl2Tracks[il][0] = 0;
+ if (fTCl2Tracks[il][1]) delete fTCl2Tracks[il][1]; fTCl2Tracks[il][1] = 0;
+ }
fNtracks = fNsingle = 0;
for (int i=6;i--;) fITSClusters[0] = 0;
fFiredChips[0] = fFiredChips[1] = 0;
fDeltTh[id] = tlet[3];
fLabels[id] = Int_t(tlet[4]);
fLabelsL2[id] = Int_t(tlet[5]);
- fUsedClusT[id] = (((ULong64_t)trSPD2)<<32) + trSPD1;
+ if (!GetMultTrackRefs()) fUsedClusT[id] = (((ULong64_t)trSPD2)<<32) + trSPD1;
//
}
fThsingle[id] = scl[0];
fPhisingle[id] = scl[1];
fLabelssingle[id] = Int_t(scl[2]);
- fUsedClusS[id] = tr;
+ if (!GetMultTrackRefs()) fUsedClusS[id] = tr;
//
}
Bool_t AliMultiplicity::FreeClustersTracklet(Int_t i, Int_t mode) const
{
// return kTRUE if the tracklet was not used by the track (on any of layers) of type mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure
- if (!fUsedClusT || mode<0 || mode>1 || i<0 || i>fNtracks) return kFALSE;
+ if (mode<0 || mode>1 || i<0 || i>fNtracks) return kFALSE;
+ if (GetMultTrackRefs()) { // new format allows multiple references
+ return !((fTCl2Tracks[0][mode] && fTCl2Tracks[0][mode]->HasReference(i)) ||
+ (fTCl2Tracks[1][mode] && fTCl2Tracks[1][mode]->HasReference(i)));
+ }
+ //
+ if (!fUsedClusT) return kFALSE;
const ULong64_t kMask0 = 0x0000ffff0000ffffLL;
const ULong64_t kMask1 = 0xffff0000ffff0000LL;
return (fUsedClusT[i]&(mode ? kMask1:kMask0)) == 0;
// set spd1 and spd2 to ID's of the tracks using the clusters of the tracklet (-1 if not used)
// Mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure tracks
// return false if the neither of clusters is used
- //
// note: stored value: [(idSAPureSPD2+1)<<16+(idTPCITS/SA_SPD2+1)]<<32 + [(idSAPureSPD1+1)<<16+(idTPCITS/SA_SPD1+1)]
- if (!fUsedClusT || mode<0 || mode>1 || i<0 || i>fNtracks) {spd1 = spd2 = -1; return kFALSE;}
- //
- spd1 = (fUsedClusT[i]&0xffffffffLL);
- spd2 = (fUsedClusT[i]>>32);
- if (mode) {
- spd1 >>= 16;
- spd2 >>= 16;
+ // Attention: new format allows references to multiple tracks, here only the 1st will be returned
+ spd1 = spd2 = -1;
+ if ( mode<0 || mode>1 || i<0 || i>fNtracks ) return kFALSE;
+ if (GetMultTrackRefs()) {
+ if (fTCl2Tracks[0][mode]) spd1 = fTCl2Tracks[0][mode]->GetReference(i,0);
+ if (fTCl2Tracks[1][mode]) spd2 = fTCl2Tracks[1][mode]->GetReference(i,0);
}
else {
- spd1 &= 0xffff;
- spd2 &= 0xffff;
+ if (!fUsedClusT) return kFALSE;
+ spd1 = (fUsedClusT[i]&0xffffffffLL);
+ spd2 = (fUsedClusT[i]>>32);
+ if (mode) { spd1 >>= 16; spd2 >>= 16;}
+ else { spd1 &= 0xffff; spd2 &= 0xffff;}
+ spd1--; // we are storing id+1
+ spd2--;
}
- spd1--; // we are storing id+1
- spd2--;
return !(spd1<0&&spd2<0);
}
+//______________________________________________________________________
+Int_t AliMultiplicity::GetTrackletTrackIDsLay(Int_t lr,Int_t i, Int_t mode, UInt_t* refs, UInt_t maxRef) const
+{
+ // fill array refs with maximum maxRef references on tracks used by the cluster of layer lr of tracklet i.
+ // return number of filled references
+ // Mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure tracks
+ //
+ int nrefs = 0;
+ if ( mode<0 || mode>1 || i<0 || i>fNtracks || lr<0||lr>1) return nrefs;
+ if (GetMultTrackRefs()) {
+ if (fTCl2Tracks[lr][mode]) nrefs = fTCl2Tracks[lr][mode]->GetReferences(i,refs, maxRef);
+ }
+ else {
+ if (!fUsedClusT || maxRef<1) return nrefs;
+ int tr = (fUsedClusT[i]&0xffffffffLL);
+ if (mode) { lr==0 ? tr >>= 16 : tr >>= 16;}
+ else { lr==0 ? tr &= 0xffff : tr &= 0xffff;}
+ refs[0] = tr--; // we are storing id+1
+ nrefs = 1;
+ }
+ return nrefs;
+}
+
+//______________________________________________________________________
+Int_t AliMultiplicity::GetSingleClusterTrackIDs(Int_t i, Int_t mode, UInt_t* refs, UInt_t maxRef) const
+{
+ // fill array refs with maximum maxRef references on tracks used by the single cluster i of layer lr
+ // return number of filled references
+ // Mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure tracks
+ //
+ int nrefs = 0;
+ if ( mode<0 || mode>1 || i<0 || i>fNtracks) return nrefs;
+ if (GetMultTrackRefs()) {
+ if (fSCl2Tracks[mode]) nrefs = fSCl2Tracks[mode]->GetReferences(i,refs, maxRef);
+ }
+ else {
+ if (!fUsedClusS || maxRef<1) return nrefs;
+ int tr = fUsedClusS[i];
+ if (mode) tr >>= 16;
+ else tr &= 0xffff;
+ refs[0] = tr--; // we are storing id+1
+ nrefs = 1;
+ }
+ return nrefs;
+}
+
//______________________________________________________________________
Bool_t AliMultiplicity::FreeSingleCluster(Int_t i, Int_t mode) const
{
// return kTRUE if the cluster was not used by the track of type mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure
- if (!fUsedClusS || mode<0 || mode>1 || i<0 || i>fNsingle) return kFALSE;
+ if (mode<0 || mode>1 || i<0 || i>fNsingle) return kFALSE;
+ if (GetMultTrackRefs()) { // new format allows multiple references
+ return !(fSCl2Tracks[mode] && fSCl2Tracks[mode]->HasReference(i));
+ }
+ if (!fUsedClusS) return kFALSE;
const UInt_t kMask0 = 0x0000ffff;
const UInt_t kMask1 = 0xffff0000;
return (fUsedClusS[i]&(mode ? kMask1:kMask0)) == 0;
// Mode: 0=TPC/ITS or ITS_SA, 1=ITS_SA_Pure tracks
// return false if the cluster is not used
//
- // note: stored value: [(idSAPure+1)<<16+(idTPCITS/SA+1)]
- if (!fUsedClusS || mode<0 || mode>1 || i<0 || i>fNsingle) {tr = -1; return kFALSE;}
- tr = fUsedClusS[i];
- if (mode) tr >>= 16;
- else tr &= 0xffff;
- return (--tr)>=0; // we are storing id+1
+ // Attention: new format allows references to multiple tracks, here only the 1st will be returned
+ tr = -1;
+ if (mode<0 || mode>1 || i<0 || i>fNsingle) return kFALSE;
+ if (GetMultTrackRefs()) { if (fSCl2Tracks[mode]) tr = fSCl2Tracks[mode]->GetReference(i,0);}
+ else {
+ if (!fUsedClusS) return kFALSE;
+ tr = fUsedClusS[i];
+ if (mode) tr >>= 16; else tr &= 0xffff;
+ tr--;
+ }
+ return tr>=0;
}
//______________________________________________________________________
void AliMultiplicity::Print(Option_t *opt) const
{
// print
- printf("N.tracklets: %4d N.singles: %4d\n",fNtracks,fNsingle);
+ printf("N.tracklets: %4d N.singles: %4d, Multiple cluster->track refs:%s\n"
+ "Used: DPhiShift: %.3e Sig^2: dPhi:%.3e dTht:%.3e NStdDev:%.2f ScaleDThtSin2T:%s\n",
+ fNtracks,fNsingle,GetMultTrackRefs() ? "ON":"OFF",
+ fDPhiShift,fDPhiWindow2,fDThetaWindow2,fNStdDev,GetScaleDThetaBySin2T() ? "ON":"OFF");
TString opts = opt; opts.ToLower();
int t0spd1=-1,t1spd1=-1,t0spd2=-1,t1spd2=-1;
//
#ifndef ALIMULTIPLICITY_H
#define ALIMULTIPLICITY_H
-#include<TObject.h>
+#include <TObject.h>
#include <TBits.h>
-#include<TMath.h>
+#include <TMath.h>
+class AliRefArray;
////////////////////////////////////////////////////////
//// Class containing multiplicity information //
class AliMultiplicity : public TObject {
public:
-
+ //
+ enum {kMultTrackRefs =BIT(14),// in new format (old is default for bwd.comp.) multiple cluster->track references are allowed
+ kScaleDThtbySin2=BIT(15) // scale Dtheta by 1/sin^2(theta). Default is DON'T scale, for bwd.comp.
+ };
AliMultiplicity(); // default constructor
// standard constructor
AliMultiplicity(Int_t ntr,Float_t *th, Float_t *ph, Float_t *dth, Float_t *dph, Int_t *labels,
virtual void Copy(TObject &obj) const;
virtual void Clear(Option_t* opt="");
virtual ~AliMultiplicity();
-// methods to access tracklet information
+ // methods to access tracklet information
+ Bool_t GetMultTrackRefs() const {return TestBit(kMultTrackRefs);}
+ void SetMultTrackRefs(Bool_t v) {SetBit(kMultTrackRefs,v);}
+ Bool_t GetScaleDThetaBySin2T() const {return TestBit(kScaleDThtbySin2);}
+ void SetScaleDThetaBySin2T(Bool_t v) {SetBit(kScaleDThtbySin2,v);}
+
+ //
Int_t GetNumberOfTracklets() const {return fNtracks;}
Double_t GetTheta(Int_t i) const {
if(i>=0 && i<fNtracks) return fTh[i];
Error("GetDeltaPhi","Invalid track number %d",i); return -9999.;
}
+ Double_t CalcDist(Int_t it) const;
+
Int_t GetLabel(Int_t i, Int_t layer) const;
void SetLabel(Int_t i, Int_t layer, Int_t label);
Int_t GetLabelSingle(Int_t i) const;
void SetLabelSingle(Int_t i, Int_t label);
-
Bool_t FreeClustersTracklet(Int_t i, Int_t mode) const;
Bool_t FreeSingleCluster(Int_t i, Int_t mode) const;
Bool_t TestFiredChipMap(UInt_t chipKey) const {return fClusterFiredChips.TestBitNumber(chipKey);}
Bool_t GetTrackletTrackIDs(Int_t i, Int_t mode, Int_t &spd1, Int_t &spd2) const;
+ Int_t GetTrackletTrackIDsLay(Int_t lr,Int_t i, Int_t mode, UInt_t* refs, UInt_t maxRef) const;
Bool_t GetSingleClusterTrackID(Int_t i, Int_t mode, Int_t &tr) const;
+ Int_t GetSingleClusterTrackIDs(Int_t i, Int_t mode, UInt_t* refs, UInt_t maxRef) const;
// array getters
Double_t* GetTheta() const {return (Double_t*)fTh;}
Int_t* GetLabels() const {return (Int_t*)fLabels;}
Int_t* GetLabels2() const {return (Int_t*)fLabelsL2;}
Int_t* GetLabelsSingle() const {return (Int_t*)fLabelssingle;}
-
- void SetTrackletData(Int_t id, const Float_t* tlet, UInt_t trSPD1, UInt_t trSPD2);
- void SetSingleClusterData(Int_t id, const Float_t* scl,UInt_t tr);
- void CompactBits();
+ void AttachTracklet2TrackRefs(AliRefArray* l1t1,AliRefArray* l1t2,AliRefArray* l2t1,AliRefArray* l2t2) {
+ fTCl2Tracks[0][0] = l1t1; fTCl2Tracks[0][1] = l1t2; fTCl2Tracks[1][0] = l2t1; fTCl2Tracks[1][1] = l2t2;
+ }
+ void AttachCluster2TrackRefs(AliRefArray* l1t1,AliRefArray* l1t2) {
+ fSCl2Tracks[0] = l1t1; fSCl2Tracks[1] = l1t2;
+ }
+ void SetTrackletData(Int_t id, const Float_t* tlet, UInt_t trSPD1=0, UInt_t trSPD2=0);
+ void SetSingleClusterData(Int_t id, const Float_t* scl,UInt_t tr=0);
+ void CompactBits();
+ //
+ void SetDPhiWindow2(Float_t v=-1) {fDPhiWindow2 = v;}
+ void SetDThetaWindow2(Float_t v=-1) {fDThetaWindow2 = v;}
+ void SetDPhiShift(Float_t v=-1) {fDPhiShift = v;}
+ void SetNStdDev(Float_t v=1) {fNStdDev = v;}
+ //
+ Float_t GetDPhiWindow2() const {return fDPhiWindow2;}
+ Float_t GetDThetaWindow2() const {return fDThetaWindow2;}
+ Float_t GetDPhiShift() const {return fDPhiShift;}
+ Float_t GetNStdDev() const {return fNStdDev;}
+
+ //
virtual void Print(Option_t *opt="") const;
protected:
void Duplicate(const AliMultiplicity &m); // used by copy ctr.
Int_t fNtracks; // Number of tracklets
- Int_t fNsingle; // Number of clusters on SPD layer 1, not associated
- // with a tracklet on SPD layer 2
+ Int_t fNsingle; // Number of clusters on SPD layer 1, not associated with a tracklet on SPD layer 2
+ //
+ Float_t fDPhiWindow2; // sigma^2 in dphi used in reco
+ Float_t fDThetaWindow2; // sigma^2 in dtheta used in reco
+ Float_t fDPhiShift; // bending shift used
+ Float_t fNStdDev; // number of standard deviations kept
+ //
Int_t *fLabels; //[fNtracks] array with labels of cluster in L1 used for tracklet
Int_t *fLabelsL2; //[fNtracks] array with labels of cluster in L2 used for tracklet
- UInt_t* fUsedClusS; //[fNsingle] id+1 of the tracks using cluster, coded as (TPC/ITS+ITS_SA)+(ITS_SA_PURE<<16)
+ UInt_t* fUsedClusS; //[fNsingle] id+1 of the tracks using cluster, coded as (TPC/ITS+ITS_SA)+(ITS_SA_PURE<<16) !!! Outphased for multiple refs
ULong64_t* fUsedClusT; //[fNtracks] id+1 of the tracks using clusters, coded as (TPC/ITS+ITS_SA)+(ITS_SA_PURE<<16) for SPD1 and SPD2 in low and high parts
+ AliRefArray *fTCl2Tracks[2][2]; // container with multiple tracklet_cluster->track references
+ AliRefArray *fSCl2Tracks[2]; // container with multiple single_cluster->track references
Double32_t *fTh; //[fNtracks] array with theta values
Double32_t *fPhi; //[fNtracks] array with phi values
Double32_t *fDeltTh; //[fNtracks] array with delta theta values
}
+inline Double_t AliMultiplicity::CalcDist(Int_t i) const
+{
+ // calculate eliptical distance. theta is the angle of cl1, dtheta = tht(cl1)-tht(cl2)
+ if (i<0 && i>=fNtracks) return -1;
+ if (fDPhiWindow2<1E-9 || fDThetaWindow2<1E-9) return -1; // not stored
+ double dphi = TMath::Abs(fDeltPhi[i]) - fDPhiShift;
+ double dtheta = fDeltTh[i];
+ if (GetScaleDThetaBySin2T()) {
+ double sinTI = TMath::Sin(fTh[i]-dtheta/2);
+ sinTI *= sinTI;
+ dtheta /= sinTI>1.e-6 ? sinTI : 1.e-6;
+ }
+ return dphi*dphi/fDPhiWindow2 + dtheta*dtheta/fDThetaWindow2;
+}
+
#endif