/** * \file ToFLevel2.cpp * \author Gianfranca DeRosa, Wolfgang Menn * * WM dec 2008: Description of "GetdEdx" changed * WM dec 2008: "GetdEdxPaddle" modified: Now includes saturation limit * PMTs higher than the saturation limit are not used for dEdx * WM apr 2009: bug found by Nicola in method "GetPaddlePlane" */ #include using namespace std; ClassImp(ToFPMT); ClassImp(ToFdEdx); ClassImp(ToFGeom); ClassImp(ToFTrkVar); ClassImp(ToFLevel2); ToFPMT::ToFPMT(){ pmt_id = 0; adc = 0.; tdc_tw = 0.; tdc = 0.; l0flag_adc = 0.; l0flag_tdc = 0.; } ToFPMT::ToFPMT(const ToFPMT &t){ pmt_id = t.pmt_id; adc = t.adc; tdc_tw = t.tdc_tw; tdc = t.tdc; } void ToFPMT::Clear(Option_t *t){ pmt_id = 0; adc = 0.; tdc_tw = 0.; tdc = 0.; } ToFTrkVar::ToFTrkVar() { trkseqno = 0; npmttdc = 0; npmtadc = 0; pmttdc = TArrayI(48); pmtadc = TArrayI(48); tdcflag = TArrayI(48); // gf: 30 Nov 2006 adcflag = TArrayI(48); // gf: 30 Nov 2006 dedx = TArrayF(48); // // memset(beta, 0, 13*sizeof(Float_t)); memset(xtofpos, 0, 3*sizeof(Float_t)); memset(ytofpos, 0, 3*sizeof(Float_t)); memset(xtr_tof, 0, 6*sizeof(Float_t)); memset(ytr_tof, 0, 6*sizeof(Float_t)); // }; void ToFTrkVar::Clear(Option_t *t) { trkseqno = 0; npmttdc = 0; npmtadc = 0; pmttdc.Reset(); pmtadc.Reset(); tdcflag.Reset(); // gf: 30 Nov 2006 adcflag.Reset(); // gf: 30 Nov 2006 dedx.Reset(); // memset(beta, 0, 13*sizeof(Float_t)); memset(xtofpos, 0, 3*sizeof(Float_t)); memset(ytofpos, 0, 3*sizeof(Float_t)); memset(xtr_tof, 0, 6*sizeof(Float_t)); memset(ytr_tof, 0, 6*sizeof(Float_t)); // }; ToFTrkVar::ToFTrkVar(const ToFTrkVar &t){ trkseqno = t.trkseqno; // npmttdc = t.npmttdc; npmtadc = t.npmtadc; (t.pmttdc).Copy(pmttdc); (t.pmtadc).Copy(pmtadc); (t.tdcflag).Copy(tdcflag); // gf: 30 Nov 2006 (t.adcflag).Copy(adcflag); // gf: 30 Nov 2006 (t.dedx).Copy(dedx); // memcpy(beta,t.beta,sizeof(beta)); memcpy(xtofpos,t.xtofpos,sizeof(xtofpos)); memcpy(ytofpos,t.ytofpos,sizeof(ytofpos)); memcpy(xtr_tof,t.xtr_tof,sizeof(xtr_tof)); memcpy(ytr_tof,t.ytr_tof,sizeof(ytr_tof)); // }; ToFLevel2::ToFLevel2() { // // PMT = new TClonesArray("ToFPMT",12); //ELENA // ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA PMT = 0; //ELENA ToFTrk = 0; //ELENA // this->Clear(); // }; void ToFLevel2::Set(){//ELENA if(!PMT)PMT = new TClonesArray("ToFPMT",12); //ELENA if(!ToFTrk)ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA }//ELENA //-------------------------------------- // // //-------------------------------------- void ToFLevel2::SetTrackArray(TClonesArray *track){//ELENA if(track && strcmp(track->GetClass()->GetName(),"ToFTrkVar")==0){ if(ToFTrk)ToFTrk->Clear("C"); ToFTrk = track; } } void ToFLevel2::Clear(Option_t *t){ // if(ToFTrk)ToFTrk->Delete(); //ELENA if(PMT)PMT->Delete(); //ELENA memset(tof_j_flag, 0, 6*sizeof(Int_t)); unpackError = 0; unpackWarning = 0; // }; void ToFLevel2::Delete(Option_t *t){ //ELENA // if(ToFTrk){ ToFTrk->Delete(); //ELENA delete ToFTrk; //ELENA } if(PMT){ PMT->Delete(); //ELENA delete PMT; //ELENA } //ELENA // }; //ELENA ToFTrkVar *ToFLevel2::GetToFTrkVar(Int_t itrk){ // if(itrk >= ntrk()){ printf(" ToFLevel2 ERROR: track related variables set %i does not exists! \n",itrk); printf(" stored track related variables = %i \n",ntrk()); return(NULL); } // if(!ToFTrk)return 0; //ELENA TClonesArray &t = *(ToFTrk); ToFTrkVar *toftrack = (ToFTrkVar*)t[itrk]; return toftrack; } /** * Retrieves the tof track matching the seqno-th tracker stored track. * */ ToFTrkVar *ToFLevel2::GetToFStoredTrack(int seqno){ if( ntrk()==0 ){ printf("ToFLevel2::GetToFStoredTrack(int) : requested tracker SeqNo %i but no ToFrimeter tracks are stored\n",seqno); return NULL; }; ToFTrkVar *c = 0; Int_t it_tof=0; do { c = GetToFTrkVar(it_tof); it_tof++; } while( c && seqno != c->trkseqno && it_tof < ntrk()); if(!c || seqno != c->trkseqno){ c = 0; if(seqno!=-1 ) printf("ToFLevel2::GetToFStoredTrack(int) : requested tracker SeqNo %i does not match ToFrimeter stored tracks\n",seqno); }; return c; } ToFPMT *ToFLevel2::GetToFPMT(Int_t ihit){ // if(ihit >= npmt()){ printf(" ToFLevel2 ERROR: pmt variables set %i does not exists! \n",ihit); printf(" stored pmt variables = %i \n",npmt()); return(NULL); } // if(!PMT)return 0; //ELENA TClonesArray &t = *(PMT); ToFPMT *tofpmt = (ToFPMT*)t[ihit]; return tofpmt; } //-------------------------------------- // // //-------------------------------------- /** * Method to get the plane ID (11 12 21 22 31 32) from the plane index (0 1 2 3 4 5) * @param Plane index (0,1,2,3,4,5). */ Int_t ToFLevel2::GetToFPlaneID(Int_t ip){ if(ip>=0 && ip<6)return 10*((int)(ip/2+1.1))+(ip%2)+1; else return -1; }; /** * Method to get the plane index (0 1 2 3 4 5) from the plane ID (11 12 21 22 31 32) * @param plane Plane ID (11, 12, 21, 22, 31, 32) */ Int_t ToFLevel2::GetToFPlaneIndex(Int_t plane_id){ if( plane_id == 11 || plane_id == 12 || plane_id == 21 || plane_id == 22 || plane_id == 31 || plane_id == 32 || false)return (Int_t)(plane_id/10)*2-1- plane_id%2; else return -1; }; /** * Method to know if a given ToF paddle was hit, that is there is a TDC signal * from both PMTs. The method uses the "tof_j_flag" variable. * @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. */ Bool_t ToFLevel2::HitPaddle(Int_t plane, Int_t paddle_id){ //<<< NEW Int_t ip = -1; if (plane>=6 ) ip = GetToFPlaneIndex(plane); else if(plane>=0 && plane < 6) ip = plane; Int_t flag=0; if(ip != -1)flag = tof_j_flag[ip] & (int)pow(2.,(double)paddle_id); if( (ip == 0 && paddle_id < 8 && flag) || (ip == 1 && paddle_id < 6 && flag) || (ip == 2 && paddle_id < 2 && flag) || (ip == 3 && paddle_id < 2 && flag) || (ip == 4 && paddle_id < 3 && flag) || (ip == 5 && paddle_id < 3 && flag) || false) return true; else return false; }; /** * Strict method to get the number of hit paddles on a ToF plane. * The method uses "HitPaddle" which checks if 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). */ Int_t ToFLevel2::GetNHitPaddles(Int_t plane){ Int_t npad=0; for(Int_t i=0; i<8; i++)npad = npad + (int)HitPaddle(plane,i); return npad; }; /** * Optional method to get the number of hit paddles on a ToF plane. * The method does NOT check if there is a signal from both PMTs, it only * checks if there is some PMT signal in a paddle * @param plane Plane ID (11, 12, 21, 22, 31, 32) or Plane index (0,1,2,3,4,5). */ Int_t ToFLevel2::GetTrueNHitPaddles(Int_t plane){ Int_t npad=0; TClonesArray* Pmt = this->PMT; int paddle[24]; memset(paddle,0, 24*sizeof(int)); for(int i=0; iGetEntries(); i++) { //loop per vedere quale TOF è colpito ToFPMT* pmthit = (ToFPMT*)Pmt->At(i); int pplane = -1; int ppaddle = -1; GetPMTPaddle(pmthit->pmt_id,pplane,ppaddle); if ( pplane == plane ) paddle[ppaddle]++; } for(int i=0;i<24;i++) if ( paddle[i]>0 ) npad++; return npad; }; //new, wm Feb 15 //wm Nov 08 //gf Apr 07 /** * Method to get the mean dEdx from a ToF layer * By definition there should be PMTs with dEdx values only in one paddle of a layer * (the paddle hitted by the track), this method looks for the hitted paddle * and gives the mean dEdx of that paddle as the output * The method was modified for the "ToF-standalone" part in february 2015 * The "adcfl" option is not very useful (an artificial dEdx is per * definition= 1 mip and not a real measurement), anyway left in the code * @param notrack Track Number * @param plane Plane index (0,1,2,3,4,5) * @param adcflag in the plane (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; ) */ Float_t ToFLevel2::GetdEdx(Int_t notrack, Int_t plane, Int_t adcfl){ ToFTrkVar *trk = GetToFTrkVar(notrack); return this->GetdEdx(trk, plane, adcfl); } //new, wm Feb 15 //wm Nov 08 //gf Apr 07 /** * Method to get the mean dEdx from a ToF layer * By definition there should be PMTs with dEdx values only in one paddle of a layer * (the paddle hitted by the track), this method looks for the hitted paddle * and gives the mean dEdx of that paddle as the output * The method was modified for the "ToF-standalone" part in february 2015 * The "adcfl" option is not very useful (an artificial dEdx is per * definition= 1 mip and not a real measurement), anyway left in the code * @param trk Pointer to TofTrkVar object * @param plane Plane index (0,1,2,3,4,5) * @param adcflag in the plane (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; ) */ Float_t ToFLevel2::GetdEdx(ToFTrkVar *trk, Int_t plane, Int_t adcfl){ Float_t dedx = 0.; Float_t PadEdx =0.; Int_t SatWarning; Int_t pad=-1; // if(!trk) return 0; //ELENA // // ToF standalone part // if ( trk->trkseqno == -1 ){ // ToFTrkVar *t_tof = trk; // Find the hitted paddle (two good TDC values) using the tof_j_flag (from tofl2com.for) Int_t Ipaddle=-1; // if tof_j_flag == 0: no paddle was hitted. Otherwise decode tof_j_flag to get the paddle if (this->tof_j_flag[plane] > 0) Ipaddle = (Int_t)log2(this->tof_j_flag[plane]) ; Ipaddle = (Int_t)log2(this->tof_j_flag[plane]) ; // Get the dEdx of this paddle using "GetdEdxPaddle" if (Ipaddle>-1) { Int_t pad = GetPaddleid(plane,Ipaddle); GetdEdxPaddle(trk, pad, adcfl, PadEdx, SatWarning); dedx = PadEdx; } // If there was no correct hitted paddle, but there was one (and only one) paddle with some // PMT entries in the PMT-class (found with "GetTrueNHitPaddles", use the dEdx of this paddle if ((Ipaddle<0) && (GetTrueNHitPaddles(plane)==1)) { // find the paddle by looping over the paddles in each layer // since GetTrueNHitPaddles==1 this is OK for (Int_t ii=0; iinpmtadc; i++){ // pmt_id = (trk->pmtadc).At(i); // GetPMTIndex(pmt_id,hh,kk); adc[kk][hh] = (trk->dedx).At(i); // }; // for (Int_t i=0; ipmt_id,hh,kk); // tdc[kk][hh] = pmt->tdc_tw; // }; // return; }; /** * Method to get the plane index (0 - 5) for the PMT_ID as input * @param pmt_id PMT_ID (0 - 47) */ Int_t ToFLevel2::GetPlaneIndex(Int_t pmt_id){ TString pmtname = GetPMTName(pmt_id); pmtname.Resize(3); if ( !strcmp(pmtname,"S11") ) return(0); if ( !strcmp(pmtname,"S12") ) return(1); if ( !strcmp(pmtname,"S21") ) return(2); if ( !strcmp(pmtname,"S22") ) return(3); if ( !strcmp(pmtname,"S31") ) return(4); if ( !strcmp(pmtname,"S32") ) return(5); return(-1); }; /** * Method to get the PMT_ID if the index (4,12) is given. We have 4 channels on * each of the 12 half-boards, this method decodes which PMT is cables to which * channel. * @param hh Channel * @param kk HalfBoard */ Int_t ToFLevel2::GetPMTid(Int_t hh, Int_t kk){ // static const short tof[4][24] = { {4, 4, 4, 4, 1, 1, 2, 2, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 2, 3, 3, 3, 3, 4}, {1, 3, 5, 7, 10, 12, 2, 4, 2, 4, 6, 8, 10, 12, 1, 5, 3, 9, 7, 9, 11, 1, 5, 9}, {2, 2, 2, 2, 1, 1, 1, 1, 4, 4, 4, 4, 4, 4, 2, 1, 2, 1, 2, 2, 2, 3, 3, 4}, {6, 8, 12, 10, 8, 6, 4, 2, 12, 10, 8, 6, 4, 2, 9, 7, 11, 11, 5, 3, 1, 3, 7, 11} }; // Int_t ind = 0; Int_t k = 0; while (k < 24){ Int_t j = 0; while (j < 2){ Int_t ch = tof[2*j][k] - 1; Int_t hb = tof[2*j + 1][k] - 1; /* tofEvent->tdc[ch][hb] */ if( ch == hh && hb == kk ){ ind = 2*k + j; break; }; j++; }; k++; }; return ind; }; /** * Method to get the PMT index if the PMT ID is given. This method is the * "reverse" of method "GetPMTid" * @param ind PMT_ID (0 - 47) * @param hb HalfBoard * @param ch Channel */ void ToFLevel2::GetPMTIndex(Int_t ind, Int_t &hb, Int_t &ch){ // static const short tof[4][24] = { {4, 4, 4, 4, 1, 1, 2, 2, 3, 3, 3, 3, 3, 3, 1, 1, 1, 1, 2, 3, 3, 3, 3, 4}, {1, 3, 5, 7, 10, 12, 2, 4, 2, 4, 6, 8, 10, 12, 1, 5, 3, 9, 7, 9, 11, 1, 5, 9}, {2, 2, 2, 2, 1, 1, 1, 1, 4, 4, 4, 4, 4, 4, 2, 1, 2, 1, 2, 2, 2, 3, 3, 4}, {6, 8, 12, 10, 8, 6, 4, 2, 12, 10, 8, 6, 4, 2, 9, 7, 11, 11, 5, 3, 1, 3, 7, 11} }; // Int_t k = 0; while (k < 24){ Int_t j = 0; while (j < 2){ /* tofEvent->tdc[ch][hb] */ if( ind == 2*k + j ){ ch = tof[2*j][k] - 1; hb = tof[2*j + 1][k] - 1; return; }; j++; }; k++; }; return; }; // wm Nov 08 revision - saturation values included /// gf Apr 07 /** * Method to get the dEdx from a given ToF paddle. * If two PMTs are good, the mean dEdx of both PMTs is taken, otherwise * just the dEdx of the "good" PMT. If both PMTs are above saturation => dEdx=1000 * @param notrack Track Number * @param Paddle index (0,1,...,23). * @param adcflag in the paddle (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; ) * @param PadEdx dEdx from a given ToF paddle * @param SatWarning 1 if the PMT ios near saturation region (adcraw ~3000) */ void ToFLevel2::GetdEdxPaddle(Int_t notrack, Int_t paddleid, Int_t adcfl, Float_t &PadEdx, Int_t &SatWarning){ ToFTrkVar *trk = GetToFTrkVar(notrack); this->GetdEdxPaddle(trk, paddleid, adcfl, PadEdx, SatWarning); }; // // wm Nov 08 revision - saturation values included /// gf Apr 07 /** * Method to get the dEdx from a given ToF paddle. * If two PMTs are good, the mean dEdx of both PMTs is taken, otherwise * just the dEdx of the "good" PMT. If both PMTs are above saturation => dEdx=1000 * @param notrack Track Number * @param Paddle index (0,1,...,23). * @param adcflag in the paddle (100<-> independent of the adcflag; !=0&&!=100 <-> at least one PMT with adcflag!=0; ) * @param PadEdx dEdx from a given ToF paddle * @param SatWarning 1 if the PMT ios near saturation region (adcraw ~3000) */ void ToFLevel2::GetdEdxPaddle(ToFTrkVar *trk, Int_t paddleid, Int_t adcfl, Float_t &PadEdx, Int_t &SatWarning){ /* Float_t PMTsat[48] = { 3162.14, 3165.48, 3153.85, 3085.73, 3089.65, 3107.64, 3097.52, 3078.37, 3130.05, 3087.07, 3112.22, 3102.92, 3080.58, 3092.55, 3087.94, 3125.03, 3094.09, 3143.16, 3125.51, 3181.27, 3092.09, 3124.98, 3069.3, 3095.53, 3097.11, 3133.53, 3114.73, 3113.01, 3091.19, 3097.99, 3033.84, 3134.98, 3081.37, 3111.04, 3066.77, 3108.17, 3133, 3111.06, 3052.52, 3140.66, 3106.33, 3094.85, 3150.85, 3118.8, 3096.24, 3118.47,3111.36, 3117.11 } ; */ // new values from Napoli dec 2008 Float_t PMTsat[48] = { 3176.35,3178.19,3167.38,3099.73,3117.00,3126.29,3111.44,3092.27, 3146.48,3094.41,3132.13,3115.37,3099.32,3110.97,3111.80,3143.14, 3106.72,3153.44,3136.00,3188.96,3104.73,3140.45,3073.18,3106.62, 3112.48,3146.92,3127.24,3136.52,3109.59,3112.89,3045.15,3147.26, 3095.92,3121.05,3083.25,3123.62,3150.92,3125.30,3067.60,3160.18, 3119.36,3108.92,3164.77,3133.64,3111.47,3131.98,3128.87,3135.56 }; for (Int_t i=0; i<48;i++) PMTsat[i] = PMTsat[i] - 5.; // safety margin PadEdx = 0.; // SatWarning = 1000; SatWarning = 0; // 0=good, increase for each bad PMT Float_t dEdx[48] = {0}; Int_t pmt_id = -1; Float_t adcraw[48]; // if(!trk) return; //ELENA // Int_t pmtleft=-1; Int_t pmtright=-1; GetPaddlePMT(paddleid, pmtleft, pmtright); adcraw[pmtleft] = 4095; adcraw[pmtright] = 4095; for (Int_t jj=0; jjpmt_id; if(pmt_id==pmtleft){ adcraw[pmtleft] = pmt->adc; } if(pmt_id==pmtright){ adcraw[pmtright] = pmt->adc; } } for (Int_t i=0; inpmtadc; i++){ if((trk->adcflag).At(i)==0 || adcfl==100){ if((trk->pmtadc).At(i) == pmtleft)dEdx[pmtleft] = (trk->dedx).At(i); if((trk->pmtadc).At(i) == pmtright)dEdx[pmtright] = (trk->dedx).At(i); }else{ if((trk->pmtadc).At(i) == pmtleft)dEdx[pmtleft] = 0.; if((trk->pmtadc).At(i) == pmtright)dEdx[pmtright] = 0.; } } // if( adcraw[pmtleft] >3000 || adcraw[pmtright] >3000)SatWarning=1; //old version // Increase SatWarning Counter for each PMT>Sat if( adcraw[pmtleft] > PMTsat[pmtleft])SatWarning++; if( adcraw[pmtright] > PMTsat[pmtright])SatWarning++; // if ADC > sat set dEdx=1000 if( adcraw[pmtleft] > PMTsat[pmtleft]) dEdx[pmtleft] = 1000.; if( adcraw[pmtright] > PMTsat[pmtright]) dEdx[pmtright] = 1000. ; // if two PMT are good, take mean dEdx, otherwise only the good dEdx if(dEdx[pmtleft]<1000 && dEdx[pmtright]<1000) PadEdx = (dEdx[pmtleft]+dEdx[pmtright])*0.5; if(dEdx[pmtleft]==1000 && dEdx[pmtright]<1000) PadEdx = dEdx[pmtright]; if(dEdx[pmtleft]<1000 && dEdx[pmtright]==1000) PadEdx = dEdx[pmtleft]; }; // gf Apr 07 /** * Method to get the PMT name (like "S11_1A") if the PMT_ID is given. * Indexes of corresponding plane, paddle and pmt are also given as output. * @param ind PMT_ID (0 - 47) * @param iplane plane index (0 - 5) * @param ipaddle paddle index (relative to the plane) * @param ipmt pmt index (0(A), 1(B)) */ TString ToFLevel2::GetPMTName(Int_t ind, Int_t &iplane, Int_t &ipaddle,Int_t &ipmt){ TString pmtname = " "; static const TString photoS[48] = { "S11_1A", "S11_1B", "S11_2A", "S11_2B", "S11_3A", "S11_3B", "S11_4A", "S11_4B", "S11_5A", "S11_5B", "S11_6A", "S11_6B", "S11_7A", "S11_7B", "S11_8A", "S11_8B", "S12_1A", "S12_1B", "S12_2A", "S12_2B", "S12_3A", "S12_3B", "S12_4A", "S12_4B", "S12_5A", "S12_5B", "S12_6A", "S12_6B", "S21_1A", "S21_1B", "S21_2A", "S21_2B", "S22_1A", "S22_1B", "S22_2A", "S22_2B", "S31_1A", "S31_1B", "S31_2A", "S31_2B", "S31_3A", "S31_3B", "S32_1A", "S32_1B", "S32_2A", "S32_2B", "S32_3A", "S32_3B" }; pmtname = photoS[ind].Data(); TString ss = pmtname(1,2); iplane = (int)(atoi(ss.Data())/10)*2-3+atoi(ss.Data())%10; ss = pmtname(4); ipaddle = atoi(ss.Data())-1 ; if( pmtname.Contains("A") )ipmt=0; if( pmtname.Contains("B") )ipmt=1; return pmtname; }; /** * Method to get the PMT name (like "S11_1A") if the PMT_ID is given * @param ind PMT_ID (0 - 47) */ TString ToFLevel2::GetPMTName(Int_t ind){ Int_t iplane = -1; Int_t ipaddle = -1; Int_t ipmt = -1; return GetPMTName(ind,iplane,ipaddle,ipmt); }; // wm jun 08 Int_t ToFLevel2::GetPaddleIdOfTrack(Float_t xtr, Float_t ytr, Int_t plane){ return GetPaddleIdOfTrack(xtr ,ytr ,plane, 0.4); } // gf Apr 07 Int_t ToFLevel2::GetPaddleIdOfTrack(Float_t xtr, Float_t ytr, Int_t plane, Float_t margin){ Double_t xt,yt,xl,xh,yl,yh; Float_t tof11_x[8] = {-17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85}; Float_t tof12_y[6] = { -13.75,-8.25,-2.75,2.75,8.25,13.75}; Float_t tof21_y[2] = { 3.75,-3.75}; Float_t tof22_x[2] = { -4.5,4.5}; Float_t tof31_x[3] = { -6.0,0.,6.0}; Float_t tof32_y[3] = { -5.0,0.0,5.0}; // S11 8 paddles 33.0 x 5.1 cm // S12 6 paddles 40.8 x 5.5 cm // S21 2 paddles 18.0 x 7.5 cm // S22 2 paddles 15.0 x 9.0 cm // S31 3 paddles 15.0 x 6.0 cm // S32 3 paddles 18.0 x 5.0 cm Int_t paddleidoftrack=-1; // //--- S11 ------ if(plane==0){ xt = xtr; yt = ytr; paddleidoftrack=-1; yl = -33.0/2. ; yh = 33.0/2. ; if ((yt>yl)&&(ytxl)&&(xtxl)&&(xtyl)&&(ytxl)&&(xtyl)&&(ytyl)&&(ytxl)&&(xtyl)&&(ytxl)&&(xtxl)&&(xtyl)&&(yt "x"-sigma. A chi2 value is * calculated, furthermore a "quality" value by adding the weights which * are finally used. If all measurements are taken, "quality" will be = 22.47. * A chi2 cut around 3-4 and a quality-cut > 20 is needed for clean beta * measurements like antiprotons etc. * The Level2 output is derived in the fortran routines using: 10.,10.,20. * @param notrack Track Number * @param cut on residual: difference between single measurement and mean * @param cut on "quality" * @param cut on chi2 */ Float_t ToFTrkVar::CalcBeta( Float_t resmax, Float_t qualitycut, Float_t chi2cut){ Float_t bxx = 100.; // ToFTrkVar *trk = this; Float_t chi2,xhelp,beta_mean; Float_t w_i[12],quality,sw,sxw,res,betachi,beta_mean_inv; Float_t b[12],tdcfl; Int_t pmt_id,pmt_plane; for (Int_t i=0; i<12; i++){ b[i] = trk->beta[i]; } //======================================================================== //--- Find out ToF layers with artificial TDC values & fill vector --- //======================================================================== Float_t w_il[6]; for (Int_t jj=0; jj<6;jj++) { w_il[jj] = 1000.; } for (Int_t i=0; inpmttdc; i++){ // pmt_id = (trk->pmttdc).At(i); pmt_plane = ToFLevel2::GetPlaneIndex(pmt_id); tdcfl = (trk->tdcflag).At(i); if (w_il[pmt_plane] != 1.) w_il[pmt_plane] = tdcfl; //tdcflag }; //======================================================================== //--- Set weights for the 12 measurements using information for top and bottom: //--- if no measurements: weight = set to very high value=> not used //--- top or bottom artificial: weight*sqrt(2) //--- top and bottom artificial: weight*sqrt(2)*sqrt(2) //======================================================================== Int_t itop[12] = {0,0,1,1,2,2,3,3,0,0,1,1}; Int_t ibot[12] = {4,5,4,5,4,5,4,5,2,3,2,3}; xhelp= 1E09; for (Int_t jj=0; jj<12;jj++) { if (jj<4) xhelp = 0.11; // S1-S3 if ((jj>3)&&(jj<8)) xhelp = 0.18; // S2-S3 if (jj>7) xhelp = 0.28; // S1-S2 if ((w_il[itop[jj]] == 1000.) && (w_il[ibot[jj]] == 1000.)) xhelp = 1E09; if ((w_il[itop[jj]] == 1) || (w_il[ibot[jj]] == 1.)) xhelp = xhelp*1.414 ; if ((w_il[itop[jj]] == 1) && (w_il[ibot[jj]] == 1.)) xhelp = xhelp*2. ; w_i[jj] = 1./xhelp; } //======================================================================== //--- Calculate mean beta for the first time ----------------------------- //--- We are using "1/beta" since its error is gaussian ------------------ //======================================================================== Int_t icount=0; sw=0.; sxw=0.; beta_mean=100.; for (Int_t jj=0; jj<12;jj++){ if ((fabs(1./b[jj])>0.1)&&(fabs(1./b[jj])<15.)) { icount= icount+1; sxw=sxw + (1./b[jj])*w_i[jj]*w_i[jj] ; sw =sw + w_i[jj]*w_i[jj] ; } } if (icount>0) beta_mean=1./(sxw/sw); beta_mean_inv = 1./beta_mean; //======================================================================== //--- Calculate beta for the second time, use residuals of the single //--- measurements to get a chi2 value //======================================================================== icount=0; sw=0.; sxw=0.; betachi = 100.; chi2 = 0.; quality=0.; for (Int_t jj=0; jj<12;jj++){ if ((fabs(1./b[jj])>0.1)&&(fabs(1./b[jj])<15.)&&(w_i[jj]>0.01)) { res = beta_mean_inv - (1./b[jj]) ; if (fabs(res*w_i[jj])0) chi2 = chi2/(icount) ; if (icount>0) betachi=1./(sxw/sw); bxx = 100.; if ((chi2 < chi2cut)&&(quality>qualitycut)) bxx = betachi; // return(bxx); }; //////////////////////////////////////////////////// //////////////////////////////////////////////////// /** * See ToFTrkVar::CalcBeta(Float_t,Float_t, Float_t). */ Float_t ToFLevel2::CalcBeta(Int_t notrack, Float_t resmax, Float_t qualitycut, Float_t chi2cut){ // cout<<" in CalcBeta "<CalcBeta(resmax,qualitycut,chi2cut); }; //////////////////////////////////////////////////// //////////////////////////////////////////////////// /** * Fills a struct cToFLevel2 with values from a ToFLevel2 object (to put data into a F77 common). */ void ToFLevel2::GetLevel2Struct(cToFLevel2 *l2) const{ for(Int_t i=0;i<6;i++) l2->tof_j_flag[i]=tof_j_flag[i]; if(ToFTrk){ //ELENA l2->ntoftrk = ToFTrk->GetEntries(); for(Int_t j=0;jntoftrk;j++){ l2->toftrkseqno[j]= ((ToFTrkVar*)ToFTrk->At(j))->trkseqno; l2->npmttdc[j]= ((ToFTrkVar*)ToFTrk->At(j))->npmttdc; for(Int_t i=0;inpmttdc[j];i++){ l2->pmttdc[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->pmttdc.At(i); l2->tdcflag[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->tdcflag.At(i); // gf: 30 Nov 2006 } for(Int_t i=0;i<13;i++) l2->beta[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->beta[i]; l2->npmtadc[j]= ((ToFTrkVar*)ToFTrk->At(j))->npmtadc; for(Int_t i=0;inpmtadc[j];i++){ l2->pmtadc[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->pmtadc.At(i); l2->adcflag[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->adcflag.At(i); // gf: 30 Nov 2006 l2->dedx[i][j] = ((ToFTrkVar*)ToFTrk->At(j))->dedx.At(i); } for(Int_t i=0;i<3;i++){ l2->xtofpos[i][j]=((ToFTrkVar*)ToFTrk->At(j))->xtofpos[i]; l2->ytofpos[i][j]=((ToFTrkVar*)ToFTrk->At(j))->ytofpos[i]; } for(Int_t i=0;i<6;i++){ l2->xtr_tof[i][j]=((ToFTrkVar*)ToFTrk->At(j))->xtr_tof[i]; l2->ytr_tof[i][j]=((ToFTrkVar*)ToFTrk->At(j))->ytr_tof[i]; } } } //ELENA if(PMT){ //ELENA l2->npmt = PMT->GetEntries(); for(Int_t j=0;jnpmt;j++){ l2->pmt_id[j] = ((ToFPMT*)PMT->At(j))->pmt_id; l2->adc[j] =((ToFPMT*)PMT->At(j))->adc; l2->tdc_tw[j] =((ToFPMT*)PMT->At(j))->tdc_tw; } } //ELENA } // // Reprocessing tool // Emiliano 08/04/07 // Int_t ToFLevel2::Process(TrkLevel2 *trk, TrigLevel2 *trg, GL_RUN *run, OrbitalInfo *orb, Bool_t force){ // // Copiare qui qualcosa di simile a calonuclei per evitare di riprocessare sempre tutto // printf("\n\n\n ERROR: NOT IMPLEMENTED ANYMORE, write Emiliano if you need this method (Emiliano.Mocchiutti@ts.infn.it) \n\n\n"); return(-1); // // // // structures to communicate with F77 // // // extern struct ToFInput tofinput_; // extern struct ToFOutput tofoutput_; // // // // DB connection // // // TString host; // TString user; // TString psw; // const char *pamdbhost=gSystem->Getenv("PAM_DBHOST"); // const char *pamdbuser=gSystem->Getenv("PAM_DBUSER"); // const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW"); // if ( !pamdbhost ) pamdbhost = ""; // if ( !pamdbuser ) pamdbuser = ""; // if ( !pamdbpsw ) pamdbpsw = ""; // if ( strcmp(pamdbhost,"") ) host = pamdbhost; // if ( strcmp(pamdbuser,"") ) user = pamdbuser; // if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw; // // // // // TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data()); // if ( !dbc->IsConnected() ) return 1; // stringstream myquery; // myquery.str(""); // myquery << "SET time_zone='+0:00';"; // dbc->Query(myquery.str().c_str()); // delete dbc->Query("SET sql_mode = 'NO_UNSIGNED_SUBTRACTION';"); // GL_PARAM *glparam = new GL_PARAM(); // glparam->Query_GL_PARAM(1,1,dbc); // parameters stored in DB in GL_PRAM table // trk->LoadField(glparam->PATH+glparam->NAME); // // // Bool_t defcal = true; // Int_t error=glparam->Query_GL_PARAM(run->RUNHEADER_TIME,201,dbc); // parameters stored in DB in GL_PRAM table // if ( error<0 ) { // return(1); // }; // printf(" Reading ToF parameter file: %s \n",(glparam->PATH+glparam->NAME).Data()); // if ( (UInt_t)glparam->TO_TIME != (UInt_t)4294967295UL ) defcal = false; // // // Int_t nlen = (Int_t)(glparam->PATH+glparam->NAME).Length(); // rdtofcal((char *)(glparam->PATH+glparam->NAME).Data(),&nlen); // // // Int_t adc[4][12]; // Int_t tdc[4][12]; // Float_t tdcc[4][12]; // // // // process tof data // // // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // adc[kk][hh] = 4095; // tdc[kk][hh] = 4095; // tdcc[kk][hh] = 4095.; // tofinput_.adc[hh][kk] = 4095; // tofinput_.tdc[hh][kk] = 4095; // }; // }; // Int_t ntrkentry = 0; // Int_t npmtentry = 0; // Int_t gg = 0; // Int_t hh = 0; // Int_t adcf[48]; // memset(adcf, 0, 48*sizeof(Int_t)); // Int_t tdcf[48]; // memset(tdcf, 0, 48*sizeof(Int_t)); // for (Int_t pm=0; pm < this->ntrk() ; pm++){ // ToFTrkVar *ttf = this->GetToFTrkVar(pm); // for ( Int_t nc=0; nc < ttf->npmttdc; nc++){ // if ( (ttf->tdcflag).At(nc) != 0 ) tdcf[(ttf->pmttdc).At(nc)] = 1; // }; // for ( Int_t nc=0; nc < ttf->npmtadc; nc++){ // if ( (ttf->adcflag).At(nc) != 0 ) adcf[(ttf->pmtadc).At(nc)] = 1; // }; // }; // // // for (Int_t pm=0; pm < this->npmt() ; pm++){ // ToFPMT *pmt = this->GetToFPMT(pm); // this->GetPMTIndex(pmt->pmt_id, gg, hh); // if ( adcf[pmt->pmt_id] == 0 ){ // tofinput_.adc[gg][hh] = (int)pmt->adc; // adc[hh][gg] = (int)pmt->adc; // }; // if ( tdcf[pmt->pmt_id] == 0 ){ // tofinput_.tdc[gg][hh] = (int)pmt->tdc; // tdc[hh][gg] = (int)pmt->tdc; // }; // tdcc[hh][gg] = (float)pmt->tdc_tw; // // Int_t pppid = this->GetPMTid(hh,gg); // // printf(" pm %i pmt_id %i pppid %i hh %i gg %i tdcc %f tdc %f adc %f \n",pm,pmt->pmt_id,pppid,hh,gg,pmt->tdc_tw,pmt->tdc,pmt->adc); // }; // // // Int_t unpackError = this->unpackError; // // // for (Int_t hh=0; hh<5;hh++){ // tofinput_.patterntrig[hh]=trg->patterntrig[hh]; // }; // // // this->Clear(); // // // Int_t pmt_id = 0; // ToFPMT *t_pmt = new ToFPMT(); // if(!(this->PMT)) this->PMT = new TClonesArray("ToFPMT",12); //ELENA // TClonesArray &tpmt = *this->PMT; // ToFTrkVar *t_tof = new ToFTrkVar(); // if(!(this->ToFTrk)) this->ToFTrk = new TClonesArray("ToFTrkVar",2); //ELENA // TClonesArray &t = *this->ToFTrk; // // // // // // Here we have calibrated data, ready to be passed to the FORTRAN routine which will extract common and track-related variables. // // // npmtentry = 0; // // // ntrkentry = 0; // // // // Calculate tracks informations from ToF alone // // // tofl2com(); // // // memcpy(this->tof_j_flag,tofoutput_.tof_j_flag,6*sizeof(Int_t)); // // // t_tof->trkseqno = -1; // // // // and now we must copy from the output structure to the level2 class: // // // t_tof->npmttdc = 0; // // // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // if ( tofoutput_.tofmask[hh][kk] != 0 ){ // pmt_id = this->GetPMTid(kk,hh); // t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc); // t_tof->tdcflag.AddAt(tofoutput_.tdcflagtof[hh][kk],t_tof->npmttdc); // gf: Jan 09/07 // t_tof->npmttdc++; // }; // }; // }; // for (Int_t kk=0; kk<13;kk++){ // t_tof->beta[kk] = tofoutput_.betatof_a[kk]; // } // // // t_tof->npmtadc = 0; // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // if ( tofoutput_.adctof_c[hh][kk] < 1000 ){ // t_tof->dedx.AddAt(tofoutput_.adctof_c[hh][kk],t_tof->npmtadc); // pmt_id = this->GetPMTid(kk,hh); // t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc); // t_tof->adcflag.AddAt(tofoutput_.adcflagtof[hh][kk],t_tof->npmtadc); // gf: Jan 09/07 // t_tof->npmtadc++; // }; // }; // }; // // // memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos)); // memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos)); // memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof)); // memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof)); // // // new(t[ntrkentry]) ToFTrkVar(*t_tof); // ntrkentry++; // t_tof->Clear(); // // // // // // // t_pmt->Clear(); // // // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // // new WM // if ( tofoutput_.tdc_c[hh][kk] < 4095 || adc[kk][hh] < 4095 || tdc[kk][hh] < 4095 ){ // // if ( tdcc[kk][hh] < 4095. || adc[kk][hh] < 4095 || tdc[kk][hh] < 4095 ){ // // // t_pmt->pmt_id = this->GetPMTid(kk,hh); // t_pmt->tdc_tw = tofoutput_.tdc_c[hh][kk]; // t_pmt->adc = (Float_t)adc[kk][hh]; // t_pmt->tdc = (Float_t)tdc[kk][hh]; // // // new(tpmt[npmtentry]) ToFPMT(*t_pmt); // npmtentry++; // t_pmt->Clear(); // }; // }; // }; // // // // Calculate track-related variables // // // if ( trk->ntrk() > 0 ){ // // // // We have at least one track // // // // // // Run over tracks // // // for(Int_t nt=0; nt < trk->ntrk(); nt++){ // // // TrkTrack *ptt = trk->GetStoredTrack(nt); // // // // Copy the alpha vector in the input structure // // // for (Int_t e = 0; e < 5 ; e++){ // tofinput_.al_pp[e] = ptt->al[e]; // }; // // // // Get tracker related variables for this track // // // toftrk(); // // // // Copy values in the class from the structure (we need to use a temporary class to store variables). // // // t_tof->npmttdc = 0; // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // if ( tofoutput_.tofmask[hh][kk] != 0 ){ // pmt_id = this->GetPMTid(kk,hh); // t_tof->pmttdc.AddAt(pmt_id,t_tof->npmttdc); // t_tof->tdcflag.AddAt(tofoutput_.tdcflag[hh][kk],t_tof->npmttdc); // gf: Jan 09/07 // t_tof->npmttdc++; // }; // }; // }; // for (Int_t kk=0; kk<13;kk++){ // t_tof->beta[kk] = tofoutput_.beta_a[kk]; // }; // // // t_tof->npmtadc = 0; // for (Int_t hh=0; hh<12;hh++){ // for (Int_t kk=0; kk<4;kk++){ // if ( tofoutput_.adc_c[hh][kk] < 1000 ){ // t_tof->dedx.AddAt(tofoutput_.adc_c[hh][kk],t_tof->npmtadc); // pmt_id = this->GetPMTid(kk,hh); // t_tof->pmtadc.AddAt(pmt_id,t_tof->npmtadc); // t_tof->adcflag.AddAt(tofoutput_.adcflag[hh][kk],t_tof->npmtadc); // gf: Jan 09/07 // t_tof->npmtadc++; // }; // }; // }; // // // memcpy(t_tof->xtofpos,tofoutput_.xtofpos,sizeof(t_tof->xtofpos)); // memcpy(t_tof->ytofpos,tofoutput_.ytofpos,sizeof(t_tof->ytofpos)); // memcpy(t_tof->xtr_tof,tofoutput_.xtr_tof,sizeof(t_tof->xtr_tof)); // memcpy(t_tof->ytr_tof,tofoutput_.ytr_tof,sizeof(t_tof->ytr_tof)); // // // // Store the tracker track number in order to be sure to have shyncronized data during analysis // // // t_tof->trkseqno = nt; // // // // create a new object for this event with track-related variables // // // new(t[ntrkentry]) ToFTrkVar(*t_tof); // ntrkentry++; // t_tof->Clear(); // // // }; // loop on all the tracks // // // this->unpackError = unpackError; // if ( defcal ){ // this->default_calib = 1; // } else { // this->default_calib = 0; // }; //}; // return(0); } bool ToFLevel2::bit(int decimal, char pos){ return( (decimal>>pos)%2 ); } bool ToFLevel2::checkPMT(TString givenpmt){ TClonesArray* Pmt = this->PMT; // printf(" ou %s entries %i \n",givenpmt.Data(),Pmt->GetEntries()); for(int i=0; iGetEntries(); i++) { ToFPMT* pmthit = (ToFPMT*)Pmt->At(i); TString pmtname = this->GetPMTName(pmthit->pmt_id); // printf(" name %s \n",pmtname.Data()); if ( !strcmp(pmtname.Data(),givenpmt.Data()) ) return true; } // printf(" PMT %s missing \n",givenpmt.Data()); return false; } bool ToFLevel2::checkPMTpatternPMThit(TrigLevel2 *trg, int &pmtpattern, int &pmtnosignal){ UInt_t *patterntrig = trg->patterntrig; pmtpattern = 0; pmtnosignal = 0; bool good = true; //S3 if ( this->bit(patterntrig[2],0) ){ pmtpattern++; if ( !this->checkPMT("S31_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],1) ){ pmtpattern++; if ( !this->checkPMT("S31_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],2) ){ pmtpattern++; if ( !this->checkPMT("S31_3A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],3) ){ pmtpattern++; if ( !this->checkPMT("S31_1B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],4) ){ pmtpattern++; if ( !this->checkPMT("S31_2B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],5) ){ pmtpattern++; if ( !this->checkPMT("S31_3B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],6) ){ pmtpattern++; if ( !this->checkPMT("S32_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],7) ){ pmtpattern++; if ( !this->checkPMT("S32_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],8) ){ pmtpattern++; if ( !this->checkPMT("S32_3A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],9) ){ pmtpattern++; if ( !this->checkPMT("S32_1B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],10) ){ pmtpattern++; if ( !this->checkPMT("S32_2B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[2],11) ){ pmtpattern++; if ( !this->checkPMT("S32_3B")){ pmtnosignal++; good = false;}} //S2 if ( this->bit(patterntrig[3],0) ){ pmtpattern++; if ( !this->checkPMT("S21_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],1) ){ pmtpattern++; if ( !this->checkPMT("S21_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],2) ){ pmtpattern++; if ( !this->checkPMT("S21_1B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],3) ){ pmtpattern++; if ( !this->checkPMT("S21_2B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],4) ){ pmtpattern++; if ( !this->checkPMT("S22_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],5) ){ pmtpattern++; if ( !this->checkPMT("S22_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],6) ){ pmtpattern++; if ( !this->checkPMT("S22_1B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[3],7) ){ pmtpattern++; if ( !this->checkPMT("S22_2B")){ pmtnosignal++; good = false;}} //S12 if ( this->bit(patterntrig[4],0) ){ pmtpattern++; if ( !this->checkPMT("S12_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],1) ){ pmtpattern++; if ( !this->checkPMT("S12_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],2) ){ pmtpattern++; if ( !this->checkPMT("S12_3A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],3) ){ pmtpattern++; if ( !this->checkPMT("S12_4A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],4) ){ pmtpattern++; if ( !this->checkPMT("S12_5A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],5) ){ pmtpattern++; if ( !this->checkPMT("S12_6A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],6) ){ pmtpattern++; if ( !this->checkPMT("S12_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],7) ){ pmtpattern++; if ( !this->checkPMT("S12_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],8) ){ pmtpattern++; if ( !this->checkPMT("S12_3A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],9) ){ pmtpattern++; if ( !this->checkPMT("S12_4B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],10) ){ pmtpattern++; if ( !this->checkPMT("S12_5B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[4],11) ){ pmtpattern++; if ( !this->checkPMT("S12_6B")){ pmtnosignal++; good = false;}} //S11 if ( this->bit(patterntrig[5],0) ){ pmtpattern++; if ( !this->checkPMT("S11_1A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],1) ){ pmtpattern++; if ( !this->checkPMT("S11_2A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],2) ){ pmtpattern++; if ( !this->checkPMT("S11_3A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],3) ){ pmtpattern++; if ( !this->checkPMT("S11_4A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],4) ){ pmtpattern++; if ( !this->checkPMT("S11_5A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],5) ){ pmtpattern++; if ( !this->checkPMT("S11_6A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],6) ){ pmtpattern++; if ( !this->checkPMT("S11_7A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],7) ){ pmtpattern++; if ( !this->checkPMT("S11_8A")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],8) ){ pmtpattern++; if ( !this->checkPMT("S11_1B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],9) ){ pmtpattern++; if ( !this->checkPMT("S11_2B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],10) ){ pmtpattern++; if ( !this->checkPMT("S11_3B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],11) ){ pmtpattern++; if ( !this->checkPMT("S11_4B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],12) ){ pmtpattern++; if ( !this->checkPMT("S11_5B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],13) ){ pmtpattern++; if ( !this->checkPMT("S11_6B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],14) ){ pmtpattern++; if ( !this->checkPMT("S11_7B")){ pmtnosignal++; good = false;}} if ( this->bit(patterntrig[5],15) ){ pmtpattern++; if ( !this->checkPMT("S11_8B")){ pmtnosignal++; good = false;}} return good; } bool ToFLevel2::checkPMTpmttrig(TrigLevel2 *trg){ // UInt_t *patterntrig = trg->patterntrig; int rS11 = 0; int rS12 = 0; int rS21 = 0; int rS22 = 0; int rS31 = 0; int rS32 = 0; // trigger configuration for the event from saved pmts TClonesArray* Pmt = this->PMT; for(int i=0; iGetEntries(); i++) { ToFPMT* pmthit = (ToFPMT*)Pmt->At(i); TString pmtname = this->GetPMTName(pmthit->pmt_id); if ( pmtname.Contains("S11") ) rS11++; if ( pmtname.Contains("S12") ) rS12++; if ( pmtname.Contains("S21") ) rS21++; if ( pmtname.Contains("S22") ) rS22++; if ( pmtname.Contains("S31") ) rS31++; if ( pmtname.Contains("S32") ) rS32++; } int rTOF1 = (rS11 + rS12) * (rS21 + rS22) * (rS31 + rS32); int rTOF2 = (rS11 * rS12) * (rS21 * rS22) * (rS31 * rS32); int rTOF3 = (rS21 + rS22) * (rS31 + rS32); int rTOF4 = (rS21 * rS22) * (rS31 * rS32); int rTOF5 = rS12 * (rS21 * rS22); int rTOF6 = (rS11 + rS12) * (rS31 + rS32); int rTOF7 = (rS11 * rS12) * (rS31 * rS32); // trigger configuration of the run bool TCTOF1 = false; bool TCTOF2 = false; bool TCTOF3 = false; bool TCTOF4 = false; bool TCTOF5 = false; bool TCTOF6 = false; bool TCTOF7 = false; if ( trg->trigconf & (1<<0) ) TCTOF1 = true; if ( trg->trigconf & (1<<1) ) TCTOF2 = true; if ( trg->trigconf & (1<<2) ) TCTOF3 = true; if ( trg->trigconf & (1<<3) ) TCTOF4 = true; if ( trg->trigconf & (1<<4) ) TCTOF5 = true; if ( trg->trigconf & (1<<5) ) TCTOF6 = true; if ( trg->trigconf & (1<<6) ) TCTOF7 = true; // do patterntrig pmts match the trigger configuration? bool pmtsconf_trigconf_match = true; if ( rTOF1 == 0 && TCTOF1 ) pmtsconf_trigconf_match = false; if ( rTOF2 == 0 && TCTOF2 ) pmtsconf_trigconf_match = false; if ( rTOF3 == 0 && TCTOF3 ) pmtsconf_trigconf_match = false; if ( rTOF4 == 0 && TCTOF4 ) pmtsconf_trigconf_match = false; if ( rTOF5 == 0 && TCTOF5 ) pmtsconf_trigconf_match = false; if ( rTOF6 == 0 && TCTOF6 ) pmtsconf_trigconf_match = false; if ( rTOF7 == 0 && TCTOF7 ) pmtsconf_trigconf_match = false; return pmtsconf_trigconf_match; } void ToFLevel2::printPMT(){ TClonesArray* Pmt = this->PMT; for(int i=0; iGetEntries(); i++) { ToFPMT* pmthit = (ToFPMT*)Pmt->At(i); TString pmtname = this->GetPMTName(pmthit->pmt_id); printf(" PMT hit: %s \n",pmtname.Data()); } } ToFdEdx::ToFdEdx() { memset(conn,0,12*sizeof(Bool_t)); memset(ts,0,12*sizeof(UInt_t)); memset(te,0,12*sizeof(UInt_t)); eDEDXpmt = new TArrayF(48); Define_PMTsat(); Clear(); } ToFdEdx::~ToFdEdx(){ Clear(); Delete(); } void ToFdEdx::Delete(Option_t *option){ if ( eDEDXpmt ){ eDEDXpmt->Set(0); if ( eDEDXpmt) delete eDEDXpmt; } } //------------------------------------------------------------------------ void ToFdEdx::CheckConnectors(UInt_t atime, GL_PARAM *glparam, TSQLServer *dbc) { for(int i=0; i<12; i++){ if(atime<=ts[i] || atime>te[i]){ Int_t error=glparam->Query_GL_PARAM(atime,210+i,dbc); // parameters stored in DB in GL_PRAM table if ( error<0 ) { conn[i]=false; ts[i]=0; te[i]=numeric_limits::max(); }; if ( !error ){ conn[i]=true; ts[i]=glparam->FROM_TIME; te[i]=glparam->TO_TIME; } if ( error>0 ){ conn[i]=false; ts[i]=glparam->TO_TIME; TSQLResult *pResult; TSQLRow *row; TString query= Form("SELECT FROM_TIME FROM GL_PARAM WHERE TYPE=%i AND FROM_TIME>=%i ORDER BY FROM_TIME ASC LIMIT 1;",210+i,atime); pResult=dbc->Query(query.Data()); if(!pResult->GetRowCount()){ te[i]=numeric_limits::max(); }else{ row=pResult->Next(); te[i]=(UInt_t)atoll(row->GetField(0)); } } // } } } //------------------------------------------------------------------------ void ToFdEdx::Clear(Option_t *option) { // // Set arrays and initialize structure // eDEDXpmt.Set(48); eDEDXpmt.Reset(-1); // Set array size and reset structure eDEDXpmt->Set(48); eDEDXpmt->Reset(-1); // Set array size and reset structure // }; //------------------------------------------------------------------------ void ToFdEdx::Print(Option_t *option) { // printf("========================================================================\n"); }; //------------------------------------------------------------------------ void ToFdEdx::Init(pamela::tof::TofEvent *tofl0) { // ToFLevel2 tf; for (Int_t gg=0; gg<4;gg++){ for (Int_t hh=0; hh<12;hh++){ // tofinput_.tdc[hh][gg]=tofEvent->tdc[gg][hh]; int mm = tf.GetPMTid(gg,hh); adc[mm]= (0xFFF & tofl0->adc[gg][hh]); // EM, exclude warning bits }; }; }; //------------------------------------------------------------------------ void ToFdEdx::Init(Int_t gg, Int_t hh, Float_t adce) { // ToFLevel2 tf; // for (Int_t gg=0; gg<4;gg++){ // for (Int_t hh=0; hh<12;hh++){ int mm = tf.GetPMTid(gg,hh); adc[mm]=adce; }; //------------------------------------------------------------------------ void ToFdEdx::Process(UInt_t atime, Float_t betamean, Float_t *xtr_tof, Float_t *ytr_tof, Int_t exitat) { bool debug = false; if ( debug ) printf(" INSIDE TOFDEDX PROCESS \n"); // the parameters should be already initialised by InitPar() // printf(" in process \n"); Clear(); // define angle: double dx = xtr_tof[1] - xtr_tof[5]; double dy = ytr_tof[0] - ytr_tof[4]; double dr = sqrt(dx*dx+dy*dy); double theta=atan(dr/76.81); // if ( xtr_tof[1] > 99. || xtr_tof[5] > 99. || ytr_tof[0] > 99. || ytr_tof[4] > 99. ) theta = 0.; for (Int_t ii=0; ii<6; ii++){ if ( xtr_tof[ii] > 99. ) xtr_tof[ii] = 0.; if ( ytr_tof[ii] > 99. ) ytr_tof[ii] = 0.; }; // if ( debug ) printf(" theta %f \n",theta); if ( debug ) printf(" xtr_tof %.1f %.1f %.1f %.1f %.1f %.1f \n",xtr_tof[0],xtr_tof[1],xtr_tof[2],xtr_tof[3],xtr_tof[4],xtr_tof[5]); if ( debug ) printf(" ytr_tof %.1f %.1f %.1f %.1f %.1f %.1f \n",ytr_tof[0],ytr_tof[1],ytr_tof[2],ytr_tof[3],ytr_tof[4],ytr_tof[5]); //--------------------- TABLE OF PERIODS WITH HV PROBLEMS ---------------------------- int Aconn=conn[0]; // PMT 0,20,22,24 int Bconn=conn[1]; // PMT 6,12,26,34 int Cconn=conn[2]; // PMT 4,14,28,32 int Dconn=conn[3]; // PMT 2,8,10,30 int Econn=conn[4]; // PMT 42,43,44,47 int Fconn=conn[5]; // PMT 7,19,23,27 int Gconn=conn[6]; // PMT 3,11,25,33 int Hconn=conn[7]; // PMT 1,9,13,21 int Iconn=conn[8]; // PMT 5,29,31,35 int Lconn=conn[9]; // PMT 37,40,45,46 int Mconn=conn[10]; // PMT 15,16,17,18 int Nconn=conn[11]; // PMT 36,38,39,41 if( false ) cout << Gconn << Iconn << Lconn <= 4095. ){ // eDEDXpmt[ii] = 0.; eDEDXpmt->AddAt(0.,ii); if ( debug ) printf(" %i adc>4095 \n",ii); continue; // EMILIANO }; if( adc[ii] >= (PMTsat[ii]-5.) && adc[ii] < 4095. ){ eDEDXpmt->AddAt(1000.,ii); if ( debug ) printf(" %i adc> pmtsat && adc<4095 \n",ii); continue; // EMILIANO }; if( adc[ii] <= 0. ) { eDEDXpmt->AddAt(1500.,ii); if ( debug ) printf(" %i adc<=0 \n",ii); continue; }; // double adcpC = f_adcPC( adc[ii] ); // - adc conversion in pC if ( exitat == 0 ){ eDEDXpmt->AddAt((Float_t)adcpC,ii); continue; } // printf(" e qua? \n"); double adccorr = adcpC*fabs(cos(theta)); if ( debug ) printf(" adccorr %f \n",adccorr); if(adccorr<=0.){ if ( debug ) printf(" %i adccorr<=0 \n",ii); // eDEDXpmt->AddAt((Float_t)adcpC,ii);//? continue; } if ( exitat == 1 ){ eDEDXpmt->AddAt((Float_t)adccorr,ii); continue; } // printf(" e quo? \n"); // int standard=0; int S115B_ok=0; int S115B_break=0; if(atime<1158720000)S115B_ok=1; else S115B_break=1; //------------------------------------------------------------------------ // printf(" e qui? \n"); //---------------------------------------------------- Z reconstruction double adcHe, adcnorm, adclin, dEdx;//, Zeta; // EM GCC4.7 adcHe=-2; adcnorm=-2; adclin=-2; dEdx=-2; // Zeta=-2;//EM GCC4.7 Double_t correction = 1.; if(Aconn==1 && (ii==0 || ii==20 || ii==22 || ii==24)){ correction = 1.675; } else if(Bconn==1 && (ii==6 || ii==12 || ii==26 || ii==34)){ correction = 2.482; } else if(Cconn==1 && (ii==4 || ii==14 || ii==28 || ii==32)){ correction = 1.464; } else if(Dconn==1 && (ii==2 || ii==8 || ii==10 || ii==30)){ correction = 1.995; } else if(Econn==1 && (ii==42 || ii==43 || ii==44 || ii==47)){ correction = 1.273; } else if(Fconn==1 && (ii==7 || ii==19 || ii==23 || ii==27)){ correction = 1.565; } else if(Mconn==1 && (ii==15 || ii==16 || ii==17 || ii==18)){ correction = 1.565; } else if(Nconn==1 && (ii==36 || ii==38 || ii==39 || ii==41)){ correction = 1.018; } else if(Hconn==1 && (ii==1 || ii==13 || ii==21 || (ii==9&&S115B_ok==1))){ correction = 1.84; } else if(S115B_break==1 && ii==9 && Hconn==1){ correction = 1.64; } else correction = 1.; if( ii==9 && S115B_break==1 ){ adcHe = f_att5B( ytr_tof[0] )/correction; } else { adcHe = Get_adc_he(ii, xtr_tof, ytr_tof)/correction; }; if(adcHe<=0){ if ( debug ) printf(" %i adcHe<=0 \n",ii); // eDEDXpmt->AddAt((Float_t)adccorr,ii); //? continue; } if ( exitat == 2 ){ if(ii==9 && S115B_break==1) eDEDXpmt->AddAt(36.*(Float_t)adccorr/adcHe,ii); else adclin = 4.*(Float_t)adccorr/adcHe; continue; } if(ii==9 && S115B_break==1) adcnorm = f_pos5B(adccorr); else adcnorm = f_pos( (parPos[ii]), adccorr); if(adcnorm<=0){ if ( debug ) printf(" %i adcnorm<=0 \n",ii); // eDEDXpmt->AddAt((Float_t)adccorr,ii);//? continue; } if ( debug ) printf(" adcnorm %f \n",adcnorm); if(ii==9 && S115B_break==1) adclin = 36.*adcnorm/adcHe; else adclin = 4.*adcnorm/adcHe; if ( debug ) printf(" adclin %f \n",adclin); if(adclin<=0){ if ( debug ) printf(" %i adclin<=0 \n",ii); // eDEDXpmt->AddAt((Float_t)adccorr,ii);//? continue; } if ( exitat == 3 ){ if(ii==9 && S115B_break==1) eDEDXpmt->AddAt((Float_t)adclin,ii); else eDEDXpmt->AddAt((Float_t)adclin,ii); continue; } // if ( betamean > 99. ){ // eDEDXpmt.AddAt((Float_t)adclin,ii); eDEDXpmt->AddAt((Float_t)adclin,ii); // printf(" AAPMT IS %i dedx is %f vector is %f \n",ii,adclin,eDEDXpmt[ii]); if ( debug ) printf(" %i betamean > 99 \n",ii); continue; }; // double dEdxHe=-2; if(ii==9 && S115B_break==1){ if( betamean <1. ) dEdxHe = f_BB5B( betamean ); else dEdxHe = 33; } else { if( betamean <1. ) dEdxHe = f_BB( (parBBneg[ii]), betamean ); else dEdxHe = parBBpos[ii]; } if ( debug ) printf(" dEdxHe %f \n",dEdxHe); if(dEdxHe<=0){ eDEDXpmt->AddAt((Float_t)adclin,ii); if ( debug ) printf(" %i dEdxHe<=0 \n",ii); continue; }; if(ii==9 && S115B_break==1) dEdx = f_desatBB5B( adclin ); else dEdx = f_desatBB((parDesatBB[ii]), adclin ); if(dEdx<=0){ eDEDXpmt->AddAt((Float_t)adclin,ii); if ( debug ) printf(" %i dEdx<=0 \n",ii); continue; }; if ( debug ) printf(" dEdx %f \n",dEdx); eDEDXpmt->AddAt((Float_t)dEdx,ii); // eDEDXpmt.AddAt((Float_t)dEdx,ii); // printf(" PMT IS %i dedx is %f vector is %f \n",ii,dEdx,eDEDXpmt[ii]); } //end loop on 48 PMT }; //------------------------------------------------------------------------ void ToFdEdx::Define_PMTsat() { Float_t sat[48] = { 3176.35,3178.19,3167.38,3099.73,3117.00,3126.29,3111.44,3092.27, 3146.48,3094.41,3132.13,3115.37,3099.32,3110.97,3111.80,3143.14, 3106.72,3153.44,3136.00,3188.96,3104.73,3140.45,3073.18,3106.62, 3112.48,3146.92,3127.24,3136.52,3109.59,3112.89,3045.15,3147.26, 3095.92,3121.05,3083.25,3123.62,3150.92,3125.30,3067.60,3160.18, 3119.36,3108.92,3164.77,3133.64,3111.47,3131.98,3128.87,3135.56 }; PMTsat.Set(48,sat); } //------------------------------------------------------------------------ void ToFdEdx::ReadParBBpos( const char *fname ) { // printf("read %s\n",fname); parBBpos.Set(48); FILE *fattin = fopen( fname , "r" ); for (int i=0; i<48; i++) { int tid=0; float tp; if(fscanf(fattin,"%d %f", &tid, &tp )!=2) break; parBBpos[i]=tp; } fclose(fattin); } //------------------------------------------------------------------------ void ToFdEdx::ReadParDesatBB( const char *fname ) { // printf("read %s\n",fname); FILE *fattin = fopen( fname , "r" ); for (int i=0; i<48; i++) { int tid=0; float tp[3]; if(fscanf(fattin,"%d %f %f %f", &tid, &tp[0], &tp[1], &tp[2] )!=4) break; parDesatBB[i].Set(3,tp); } fclose(fattin); } //------------------------------------------------------------------------ void ToFdEdx::ReadParBBneg( const char *fname ) { // printf("read %s\n",fname); FILE *fattin = fopen( fname , "r" ); for (int i=0; i<48; i++) { int tid=0; float tp[3]; if(fscanf(fattin,"%d %f %f %f", &tid, &tp[0], &tp[1], &tp[2] )!=4) break; parBBneg[i].Set(3,tp); } fclose(fattin); } //------------------------------------------------------------------------ void ToFdEdx::ReadParPos( const char *fname ) { // printf("read %s\n",fname); FILE *fattin = fopen( fname , "r" ); for (int i=0; i<48; i++) { int tid=0; float tp[4]; if(fscanf(fattin,"%d %f %f %f %f", &tid, &tp[0], &tp[1], &tp[2], &tp[3])!=5) break; parPos[i].Set(4,tp); } fclose(fattin); } //------------------------------------------------------------------------ void ToFdEdx::ReadParAtt( const char *fname ) { // printf("read %s\n",fname); FILE *fattin = fopen( fname , "r" ); for (int i=0; i<48; i++) { int tid=0; float tp[6]; if(fscanf(fattin,"%d %f %f %f %f %f %f", &tid, &tp[0], &tp[1], &tp[2], &tp[3], &tp[4], &tp[5] )!=7) break; parAtt[i].Set(6,tp); } fclose(fattin); } double ToFdEdx::f_att( TArrayF &p, float x ) { return p[0] + p[1]*x + p[2]*x*x + p[3]*x*x*x + p[4]*x*x*x*x + p[5]*x*x*x*x*x; } //------------------------------------------------------------------------ double ToFdEdx::f_att5B( float x ) { return 101.9409 + 6.643781*x + 0.2765518*x*x + 0.004617647*x*x*x + 0.0006195132*x*x*x*x + 0.00002813734*x*x*x*x*x; } double ToFdEdx::f_pos( TArrayF &p, float x ) { return p[0] + p[1]*x + p[2]*x*x + p[3]*x*x*x; } double ToFdEdx::f_pos5B( float x ) { return 15.45132 + 0.8369721*x + 0.0005*x*x; } double ToFdEdx::f_adcPC( float x ) { return 28.12+0.6312*x-5.647e-05*x*x+3.064e-08*x*x*x; } float ToFdEdx::Get_adc_he( int id, float pl_x[6], float pl_y[6]) { // // input: id - pmt [0:47} // pl_x - coord x of the tof plane // pl_y - coord y adc_he = 0; if( eGeom.GetXY(id)==1 ) adc_he = f_att( (parAtt[id]), pl_x[eGeom.GetPlane(id)] ); if( eGeom.GetXY(id)==2 ) adc_he = f_att( (parAtt[id]), pl_y[eGeom.GetPlane(id)] ); return adc_he; } //------------------------------------------------------------------------ double ToFdEdx::f_BB( TArrayF &p, float x ) { return p[0]/(x*x)*(log(x*x/(1-x*x)) - p[1]*x*x - p[2]); } //------------------------------------------------------------------------ double ToFdEdx::f_BB5B( float x ) { return 0.165797/(x*x)*(log(x*x/(1-x*x)) + 140.481*x*x + 52.9258); } //------------------------------------------------------------------------ double ToFdEdx::f_desatBB( TArrayF &p, float x ) { return p[0] + p[1]*x + p[2]*x*x; } //------------------------------------------------------------------------ double ToFdEdx::f_desatBB5B( float x ) { return -2.4 + 0.75*x + 0.009*x*x; }