BRAT 2.4.5
Class index
Full class index
brahmlib
BRAHMS
ROOT page

BrC1PidModule


class description - source file - inheritance tree

class BrC1PidModule : public BrModule


    protected:
virtual void CheckGeometry(BrChkvRdo* allHits, const BrVector3D& xy) virtual BrC1Pid* CreatePid(Float_t blobE, Float_t confL, Int_t trackId) virtual void EstimateBackground(BrChkvRdo* allHits, BrC1PidModule::BLOB& blob) virtual Float_t EstimateBlobRadius(Float_t p, Float_t beta) virtual void FindTubesHitByBlob(BrChkvRdo* allHits, const BrVector3D& xy, Float_t blobRadius, BrC1PidModule::BLOB& blob) virtual Bool_t IsTubeHitByBlob(Int_t tubeNo, const BrVector3D& xy, Float_t blobRadius) virtual Bool_t Overlap(BrC1PidModule::BLOB& blob1, BrC1PidModule::BLOB& blob2) virtual BrVector3D PointToC1BackPlane(BrDetectorTrack* t2track) virtual Int_t SeparateBlobs(BrC1PidModule::BLOB& blob1, BrC1PidModule::BLOB& blob2) virtual void Snapshot(BrChkvRdo* allHits, const BrVector3D& xy, Float_t p = 99) public:
BrC1PidModule BrC1PidModule() BrC1PidModule BrC1PidModule(const Char_t* name, const Char_t* title) BrC1PidModule BrC1PidModule(BrC1PidModule&) virtual void ~BrC1PidModule() virtual void Begin() static TClass* Class() virtual void DefineHistograms() virtual void End() virtual void Event(BrEventNode* inNode, BrEventNode* outNode) virtual void Finish() virtual void Init() virtual TClass* IsA() const virtual void Print(Option_t* option = "B") const void SetBackPlaneOffset(Float_t x = -0.307, Float_t y = -0.645) void SetBetaCut(Float_t b = 0.98) void SetBlobRadius(Float_t r = 3.94) void SetCheckGeometry(Bool_t b = kFALSE) void SetEstimateBackground(Bool_t b = kFALSE) void SetMaxNSnapshots(Int_t n = 0) void SetRefractiveIndex(Float_t n = 1.00138) void SetUseH1Beta(Bool_t b = kTRUE) virtual void ShowMembers(TMemberInspector& insp, char* parent) virtual void Streamer(TBuffer& b) void StreamerNVirtual(TBuffer& b)

Data Members

    private:
BrDetectorVolume* fC1Volume BrDetectorVolume* fT2Volume BrChkvParameters* fC1Params BrPlane3D fC1BackPlane BrVector3D fBackPlaneOffset Int_t fMaxNSnapshots max number of snapshots Int_t fNSnapshots number of snapshots taken Float_t fBlobRadius radius of light blob Bool_t fUseH1Beta to calculate blob radius Float_t fBetaCut Float_t fRefractiveIndex Float_t fPionThreshold pion momentum threshold Bool_t fCheckGeometry if on, geo check hist is created Bool_t fEstimateBackground if on, bg est. hist is created TDirectory* fHistDir TH1F* fBlobEnergy TH2F* fBlobEnergyVsP TProfile* fBlobEnergyVsPProfile TH2F* fTubeToHit TH1F* fTubeToHitX TH1F* fTubeToHitY TH1F* fBackground TH1F* fBackgroundToSignal TH1F* fNTubesInBlob TH2F* fBlobEnergyVsNTubes TH2F* fIntWithBackPlane TH1F* fOverlaps TH1F* fOverlapEnergy TH1F* fOverlapEnergyR TH1F* fConflicts TH1F* fConfLevel TH1F* fStatus TH1F* fUsedBlobRadius TList* fSnapshots TF1* fPionEVsP TF1* fKaonEVsP

Class Description

 C1 pid module - makes particle identification from C1 data.

 Finds Ffs track table and C1Rdo object in inNode. Matches tracks
 and hits in the sensitive plane of C1 and returns (called C1Pid)
 data table with BrC1Pid objects to the outNode.

 The track is pointed to the backplane of C1. The blob (of cherenkov
 light) has an overlap with a number of photo tubes (0<n<10). The
 calibrated adc signal from these tubes are added - this is called
 the blob energy.

 If two blobs (corresponding to two tracks) have overlapping tubes
 they should be separated - this method is not implemented yet.


BrC1PidModule()
 Default constructor. DO NOT USE

BrC1PidModule(const Char_t* name, const Char_t* title) : BrModule(name, title)
 Named Constructor

void DefineHistograms()
 Define histograms. They are:

void Init()
 Job-level initialisation

void Begin()
 Run-level initialisation

void Event(BrEventNode* inNode, BrEventNode* outNode)
 Requires:  ChkvRdo C1 and FfsTracks objects in inNode or outNode
 Generates: BrDataTable with BrC1Pid objects

void FindTubesHitByBlob(BrChkvRdo* allHits, const BrVector3D& xy, Float_t r, BLOB& blob)

Bool_t IsTubeHitByBlob(Int_t tubeNo, const BrVector3D& xy, Float_t r)
 finds the blob circle (center = (cX,cY), radius = fBlobRadius)
 and the tube square of tube with tubeNo. Is there and overlap?
 returns kTRUE is yes, kFALSE if no.

Bool_t Overlap(BLOB& blob1, BLOB& blob2)
 check if there are overlapping tubes in the two arrays.

Int_t SeparateBlobs(BLOB& blob1, BLOB& blob2)
 For now this method always returns -1!!!!!

Float_t EstimateBlobRadius(Float_t p, Float_t beta)
 Estimate the blob radius from the momentum (asume pion)
 r = L * tan(theta)

void CheckGeometry(BrChkvRdo* allHits, const BrVector3D& xy)
fills "dist to hit" histogram - used for geometry alignment

void EstimateBackground(BrChkvRdo* allHits, BLOB& blob)
 find n tubes that are not in blobHits
 n = number of hits in blobHits
 sum the adc's

void Snapshot(BrChkvRdo* allHits, const BrVector3D& xy, Float_t p=99)

PointToC1BackPlane(BrDetectorTrack* t2track)
 Returns the position (in local C1 coordinates) of the
 t2 track intersection with C1 backplane.
 If outside backplane returns -999,-999,-999

CreatePid(Float_t e, Float_t cl, Int_t trackId)

void End()
 Run-level finalisation

void Finish()
 Job-level finalisation

void Print(Option_t* option) const
 Print module information
 See BrModule::Print for options.
 In addition this module defines the Option:
 <fill in here>



Inline Functions


                 void SetBackPlaneOffset(Float_t x = -0.307, Float_t y = -0.645)
                 void SetBlobRadius(Float_t r = 3.94)
                 void SetRefractiveIndex(Float_t n = 1.00138)
                 void SetUseH1Beta(Bool_t b = kTRUE)
                 void SetBetaCut(Float_t b = 0.98)
                 void SetMaxNSnapshots(Int_t n = 0)
                 void SetCheckGeometry(Bool_t b = kFALSE)
                 void SetEstimateBackground(Bool_t b = kFALSE)
              TClass* Class()
              TClass* IsA() const
                 void ShowMembers(TMemberInspector& insp, char* parent)
                 void Streamer(TBuffer& b)
                 void StreamerNVirtual(TBuffer& b)
        BrC1PidModule BrC1PidModule(BrC1PidModule&)
                 void ~BrC1PidModule()

This page automatically generated by script docBrat by Christian Holm

Copyright ; 2002 BRAHMS Collaboration <brahmlib@rcf.rhic.bnl.gov>
Last Update on 2002/06/10 13:29:50 $ by ekman $

Validate HTML
Validate CSS