void ToFLevel2Ex(){ /* * * This script shows how use ToFLevel2 output without using PamLevel2 software * Notice that here no selection on image tracks is applied * All methods and all variables usage are shown * */ /* Some plot just as example */ // Beta versus Rigidity TH2F* BetaRig= new TH2F("BetaRig","Beta versus Rigidity",1000,0.1,10.,520,-1.3,1.3); // dEdx == adc_c measurement TH1F* ADChisto= new TH1F("ADChisto","ADChisto",10,0.,20.); TH1F* TestHisto = new TH1F("TestHisto","TestHisto",500,0.,5000.); gStyle->SetOptStat(111111); TrkLevel2* trk_event = new TrkLevel2(); ToFLevel2* tof_event = new ToFLevel2(); TChain *T = new TChain("Tracker"); TChain *C = new TChain("ToF"); T->Add("./8212.Level2.root"); C->Add("./8212.Level2.root"); T->SetBranchAddress("TrkLevel2", &trk_event); // set tracker event adrss C->SetBranchAddress("ToFLevel2", &tof_event); // set ToF event adrss T->AddFriend(C); // make the trees friends Int_t nevent = T->GetEntries(); cout << " Events in Tracker Tree " << nevent << endl; Int_t ntofevent = C->GetEntries(); cout << " Events in ToF Tree " << ntofevent << endl; int j = 0; int k = 0; // for (Int_t ii=0; iiGetEntry(ii); /* * infos in ToFLevel2 class * */ // j flag: it's a decoded value (array wioth 6 elements); // for each layer, it summed up if there // are more than one hit per layer // tof11_j = tof11_j + 2**(i-1) for(Int_t i=0; i<6; i++) Int_t tof11_j = tof_event->tof_j_flag[i]; /* * if you need informations from ToF and Tracker, * you have to read both the detector classes */ /* ToFTrkVar is the track related class */ for (Int_t it=0; itntrk() ; it++){ ToFTrkVar *tof = tof_event->GetToFTrkVar(it); /* trkseqno = index to get ToF alone or ToF+Tracker informations */ /* * trkseqno ==-1 -> ToF alone informations * trkseqno !=-1 correspond to the Tracker track index */ // cout << " Index " << tof->trkseqno << endl; /* Get Tracker information if it exist */ Int_t trkgood = 0; if(tof->trkseqno != -1){ TrkTrack *trk = trk_event->GetStoredTrack(tof->trkseqno); /* Select "good" Tracker tracks" */ if(trk->chi2 > 0 && trk->chi2 < 100 && trk->GetNX() >= 4 && trk->GetNY() >= 3 && trk->xgood[0] && (trk->xgood[5] || trk->xgood[4]) && trk->al[4] != 0.) { trkgood = 1; } } // Beta measurement based on ToF alone /* * beta is an array with 13 elements; * 12 correspond to the several planes combinations as described below; * the last one is a mean value */ /* * * S11 - S31 = tof->beta[0] * S11 - S32 = tof->beta[1] * S12 - S31 = tof->beta[2] * S12 - S32 = tof->beta[3] * S21 - S31 = tof->beta[4] * S21 - S32 = tof->beta[5] * S22 - S31 = tof->beta[6] * S22 - S32 = tof->beta[7] * S11 - S21 = tof->beta[8] * S11 - S22 = tof->beta[9] * S12 - S21 = tof->beta[10] * S12 - S22 = tof->beta[11] * */ Int_t ll = 0; Int_t mm = 0; /* * to get the same information for tracks from Tracker set: * if(tof->trkseqno==0) for the first tracker track; * if(tof->trkseqno==1) for the second tracker track; */ if(tof->trkseqno == -1){ for (Int_t k=0; k<13 ; k++){ Float_t Beta = tof->beta[k]; // cout << "Beta " << Beta << endl; } } /* * read the pmtid for the tdc used in evaluate beta * for the first tracker track */ Float_t tdcflagm[4][12]; for (Int_t j=0; j<4; j++){ for (Int_t k=0; k<12; k++){ tdcflagm[j][k]=0; } } if(tof->trkseqno == 0){ // npmttdc = number of pmt used for beta evaluation for (Int_t ipmtt=0; ipmttnpmttdc; ipmtt++){ // you can read the pmtid for the tdc used in evaluate beta Int_t pmttdc = tof->pmttdc[ipmtt]; // printf("pmttdc %d\n",pmttdc); Int_t tdcflag = tof->tdcflag[ipmtt]; if(tdcflag!=0) { // cout << " TDC " << tof->tdcflag[ipmtt] << " " << tdcflag << " " << pmttdc << endl; ll=0; mm=0; tof_event->GetPMTIndex(pmttdc,ll,mm); cout<<"pmttdc tdcflag "<GetPlaneIndex(pmttdc); } } ll = 0; mm = 0; /* * dEdx for each PMT with signal; * this value coprresponds to the signal from a single pmt */ Float_t adc2[4][12]; if(tof->trkseqno == 0){ // npmtadc = number of pmt used for dEdx evaluation for (Int_t ipmta=0; ipmtanpmtadc; ipmta++){ // dedx = adc_c Float_t dEdx = tof->dedx[ipmta]; // printf("dEdx %f %f\n",tof->dedx[ipmta],dEdx); // you can read the pmtid for the adc used in evaluate dedx=adc_c Int_t pmtadc = tof->pmtadc[ipmta]; // printf("pmtadc %d\n",tof->pmtadc[ipmta]); Int_t adcflag = tof->adcflag[ipmta]; // cout << " ADC " << adcflag << " " << pmtadc << " " << dEdx << endl; if(ipmta == 1) { ADChisto->Fill(tof->dedx[ipmta]); } /* to get pmt index in terms of matrix index from pmtid */ tof_event->GetPMTIndex(pmtadc,ll,mm); adc2[mm][ll] = tof->dedx[ipmta]; /* to get pmt index from matrix index */ tof_event->GetPMTid(ll,mm); Int_t ind = tof_event->GetPMTid(ll,mm); /* to get pmt name from pmtid */ tof_event->GetPMTName(pmtadc); // cout << tof_event->GetPMTName(pmtadc) << endl; } } // Position measurement based on ToF hits /* * to get the same information for tracks from Tracker set: * if(tof->trkseqno==0) for the first tracker track; * if(tof->trkseqno==1) for the second tracker track; */ if(tof->trkseqno==-1){ for (Int_t j=0; j<3; j++){ Float_t xtofpos = tof->xtofpos[j]; Float_t ytofpos = tof->ytofpos[j]; // printf("xtofpos %f %f %d\n",xtofpos,tof->xtofpos[j],ii); // printf("ytofpos %f %f %d\n",ytofpos,tof->ytofpos[j],ii); } } // Beta versus Rigidity /* * for the first tracker track tof->trkseqno == 0 * for the second tracker track tof->trkseqno == 1 * and so on... */ if(tof->trkseqno == 0 && trkgood==1){ BetaRig->Fill( trk->GetRigidity(), tof->beta[4]); } /* * Method to get the dEdx on a ToF plane. * It is evaluated as sum of all adc signals on the plane * after corrections; its usage is not advised; * This method will be modified in the next ToFLevel2 release */ if(tof->trkseqno == 0 && trkgood==1){ for (Int_t plane=0; plane<13; plane++){ Float_t MethoddEdx=tof_event->GetdEdx(it,plane); // it is the index for // tof_event->ntrk() } } // fill adc and tdc using getmatrix Float_t adc[4][12]; Float_t tdc[4][12]; tof_event->GetMatrix(it,adc,tdc); // it is the index for tof_event->ntrk$ for (Int_t j=0; j<4; j++){ for (Int_t k=0; k<12; k++){ Int_t ind = tof_event->GetPMTid(j,k); if (tdc[j][k] < 4095) cout<npmt(); Float_t adct[48]; Float_t tdct[48]; /* * ToFPMT is the PMT related class */ for (Int_t ipmt=0; ipmtnpmt() ; ipmt++){ /* tof_event->npmt() is the number of PMT with at least a TDC signal */ ToFPMT *tofpmt = tof_event->GetToFPMT(ipmt); // for each PMT with a signal, tofpmt->adc gives the raw adc adct[ipmt]=tofpmt->adc; // for each PMT with a signal, tofpmt->tdc_tw gives the tdc after Time-walk correction tdct[ipmt]=tofpmt->tdc_tw; TestHisto->Fill(adct[2]); } /* * Some useful Methods (Examples): */ Int_t nplane = 6; Int_t npaddle = 8; for (Int_t idpl=0; idpl< nplane; idpl++){ /* * Method to get the number of hit paddles on a ToF plane. * @param plane Plane ID (11, 12, 21, 22, 31, 32) or Plane index (0,1,2,3,4,5). * it is evaluated for each event in the run */ Int_t NHit = tof_event->GetNHitPaddles(idpl); // printf("Number of Hit, idPlane %d %d %d %d\n",NHit,idpl,ii,it); /** * Method to get the plane ID (11 12 21 22 31 32) from the plane index (0 1 2 3 4 5) */ Int_t ToFPlaneID=tof_event->GetToFPlaneID(idpl); // printf("ToFPlaneID, idPlane %d %d\n",ToFPlaneID,idpl); } /* * Method to know if a given ToF paddle was hit, that is there is a TDC signal * from both PMTs * @param plane Plane ID (11, 12, 21, 22, 31, 32) or Plane index (0,1,2,3,4,5). * @param paddle_id Paddle ID. * @return 1 if the paddle was hit. */ for (Int_t idpad=0; idpad< npaddle; idpad++){ Int_t HitPad = tof_event->HitPaddle(idpl,idpad); // printf("Hit Pad %d %d %d\n",HitPad,idpl,idpad); } /** * Method to get the plane index (0 1 2 3 4 5) from the plane ID (11 12 21 22 31 32) */ Int_t plane_id=11; Int_t ToFPlaneIndex = tof_event->GetToFPlaneIndex(plane_id); // printf("Plane Index %d %d\n",ToFPlaneIndex,plane_id); } // end loop over events TCanvas *BetaRigCanvas = new TCanvas("BetaRigCanvas","BetaRigCanvas", 700, 800); BetaRigCanvas->cd(); BetaRig->Draw(); TCanvas *SCanvas = new TCanvas("SCanvas","SCanvas", 700, 800); SCanvas->cd(); ADChisto->Draw(); TCanvas *TestCanvas = new TCanvas("TestCanvas","TestCanvas", 700, 800); TestCanvas->cd(); TestHisto->Draw(); }