X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=TRD%2FAliTRDCalibraFillHisto.cxx;h=558e0b0e37720cc746e32b4514b921e1dc5c6eb5;hp=bc237a0dade5efb4370be33e9b1b6f04dd22d74e;hb=2c45ad6f9fa1914ed81971a7747f73704f8f9496;hpb=db10d95af7a79d6b9a27a7a38b63186c3c0f549d diff --git a/TRD/AliTRDCalibraFillHisto.cxx b/TRD/AliTRDCalibraFillHisto.cxx index bc237a0dade..558e0b0e377 100644 --- a/TRD/AliTRDCalibraFillHisto.cxx +++ b/TRD/AliTRDCalibraFillHisto.cxx @@ -53,20 +53,21 @@ #include "AliLog.h" +#include "AliESDtrack.h" #include "AliTRDCalibraFillHisto.h" #include "AliTRDCalibraMode.h" #include "AliTRDCalibraVector.h" #include "AliTRDCalibraVdriftLinearFit.h" +#include "AliTRDCalibraExbAltFit.h" #include "AliTRDcalibDB.h" #include "AliTRDCommonParam.h" #include "AliTRDpadPlane.h" #include "AliTRDcluster.h" -#include "AliTRDtrack.h" #include "AliTRDtrackV1.h" -#include "AliTRDrawStreamBase.h" #include "AliRawReader.h" #include "AliRawReaderDate.h" #include "AliTRDgeometry.h" +#include "AliTRDfeeParam.h" #include "./Cal/AliTRDCalROC.h" #include "./Cal/AliTRDCalPad.h" #include "./Cal/AliTRDCalDet.h" @@ -140,6 +141,7 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto() ,fVector2d(kFALSE) ,fLinearFitterOn(kFALSE) ,fLinearFitterDebugOn(kFALSE) + ,fExbAltFitOn(kFALSE) ,fRelativeScale(0) ,fThresholdClusterPRF2(15.0) ,fLimitChargeIntegration(kFALSE) @@ -147,12 +149,21 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto() ,fNormalizeNbOfCluster(kFALSE) ,fMaxCluster(0) ,fNbMaxCluster(0) + ,fCutWithVdriftCalib(kFALSE) + ,fMinNbTRDtracklets(0) + ,fMinTRDMomentum(0.0) + ,fFirstRunGain(0) ,fVersionGainUsed(0) ,fSubVersionGainUsed(0) + ,fFirstRunGainLocal(0) ,fVersionGainLocalUsed(0) ,fSubVersionGainLocalUsed(0) + ,fFirstRunVdrift(0) ,fVersionVdriftUsed(0) ,fSubVersionVdriftUsed(0) + ,fFirstRunExB(0) + ,fVersionExBUsed(0) + ,fSubVersionExBUsed(0) ,fCalibraMode(new AliTRDCalibraMode()) ,fDebugStreamer(0) ,fDebugLevel(0) @@ -186,6 +197,7 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto() ,fCH2d(0x0) ,fLinearFitterArray(540) ,fLinearVdriftFit(0x0) + ,fExbAltFit(0x0) ,fCalDetGain(0x0) ,fCalROCGain(0x0) { @@ -219,6 +231,7 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c) ,fVector2d(c.fVector2d) ,fLinearFitterOn(c.fLinearFitterOn) ,fLinearFitterDebugOn(c.fLinearFitterDebugOn) + ,fExbAltFitOn(c.fExbAltFitOn) ,fRelativeScale(c.fRelativeScale) ,fThresholdClusterPRF2(c.fThresholdClusterPRF2) ,fLimitChargeIntegration(c.fLimitChargeIntegration) @@ -226,12 +239,21 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c) ,fNormalizeNbOfCluster(c.fNormalizeNbOfCluster) ,fMaxCluster(c.fMaxCluster) ,fNbMaxCluster(c.fNbMaxCluster) + ,fCutWithVdriftCalib(c.fCutWithVdriftCalib) + ,fMinNbTRDtracklets(c.fMinNbTRDtracklets) + ,fMinTRDMomentum(c.fMinTRDMomentum) + ,fFirstRunGain(c.fFirstRunGain) ,fVersionGainUsed(c.fVersionGainUsed) ,fSubVersionGainUsed(c.fSubVersionGainUsed) + ,fFirstRunGainLocal(c.fFirstRunGainLocal) ,fVersionGainLocalUsed(c.fVersionGainLocalUsed) ,fSubVersionGainLocalUsed(c.fSubVersionGainLocalUsed) + ,fFirstRunVdrift(c.fFirstRunVdrift) ,fVersionVdriftUsed(c.fVersionVdriftUsed) ,fSubVersionVdriftUsed(c.fSubVersionVdriftUsed) + ,fFirstRunExB(c.fFirstRunExB) + ,fVersionExBUsed(c.fVersionExBUsed) + ,fSubVersionExBUsed(c.fSubVersionExBUsed) ,fCalibraMode(0x0) ,fDebugStreamer(0) ,fDebugLevel(c.fDebugLevel) @@ -265,6 +287,7 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c) ,fCH2d(0x0) ,fLinearFitterArray(540) ,fLinearVdriftFit(0x0) + ,fExbAltFit(0x0) ,fCalDetGain(0x0) ,fCalROCGain(0x0) { @@ -288,6 +311,9 @@ AliTRDCalibraFillHisto::AliTRDCalibraFillHisto(const AliTRDCalibraFillHisto &c) if(c.fLinearVdriftFit){ fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit(*c.fLinearVdriftFit); } + if(c.fExbAltFit){ + fExbAltFit = new AliTRDCalibraExbAltFit(*c.fExbAltFit); + } if(c.fCalDetGain) fCalDetGain = new AliTRDCalDet(*c.fCalDetGain); if(c.fCalROCGain) fCalROCGain = new AliTRDCalROC(*c.fCalROCGain); @@ -331,6 +357,7 @@ AliTRDCalibraFillHisto::~AliTRDCalibraFillHisto() if(f) { delete f;} } if(fLinearVdriftFit) delete fLinearVdriftFit; + if(fExbAltFit) delete fExbAltFit; if (fGeo) { delete fGeo; } @@ -502,6 +529,17 @@ Bool_t AliTRDCalibraFillHisto::Init2Dhistos(Int_t nboftimebin) } } fLinearVdriftFit = new AliTRDCalibraVdriftLinearFit(); + TString nameee("Ver"); + nameee += fVersionExBUsed; + nameee += "Subver"; + nameee += fSubVersionExBUsed; + nameee += "FirstRun"; + nameee += fFirstRunExB; + nameee += "Nz"; + fLinearVdriftFit->SetNameCalibUsed(nameee); + } + if(fExbAltFitOn){ + fExbAltFit = new AliTRDCalibraExbAltFit(); } if (fPRF2dOn) { @@ -527,7 +565,7 @@ Bool_t AliTRDCalibraFillHisto::InitCalDet() // DB Setting // Get cal - AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberGainFactor",AliCDBManager::Instance()->GetRun(),fVersionGainUsed,fSubVersionGainUsed); + AliCDBEntry *entry = AliCDBManager::Instance()->Get("TRD/Calib/ChamberGainFactor",fFirstRunGain,fVersionGainUsed,fSubVersionGainUsed); if(!entry) { AliError("No gain det calibration entry found"); return kFALSE; @@ -550,6 +588,8 @@ Bool_t AliTRDCalibraFillHisto::InitCalDet() name += fVersionGainUsed; name += "Subver"; name += fSubVersionGainUsed; + name += "FirstRun"; + name += fFirstRunGain; name += "Nz"; name += fCalibraMode->GetNz(0); name += "Nrphi"; @@ -562,6 +602,8 @@ Bool_t AliTRDCalibraFillHisto::InitCalDet() namee += fVersionVdriftUsed; namee += "Subver"; namee += fSubVersionVdriftUsed; + namee += "FirstRun"; + namee += fFirstRunVdrift; namee += "Nz"; namee += fCalibraMode->GetNz(1); namee += "Nrphi"; @@ -569,6 +611,19 @@ Bool_t AliTRDCalibraFillHisto::InitCalDet() fPH2d->SetTitle(namee); + // title AliTRDCalibraVdriftLinearFit + TString nameee("Ver"); + nameee += fVersionExBUsed; + nameee += "Subver"; + nameee += fSubVersionExBUsed; + nameee += "FirstRun"; + nameee += fFirstRunExB; + nameee += "Nz"; + + + fLinearVdriftFit->SetNameCalibUsed(nameee); + + return kTRUE; @@ -611,139 +666,7 @@ Bool_t AliTRDCalibraFillHisto::InitCalPad(Int_t detector) } //____________Offline tracking in the AliTRDtracker____________________________ -Bool_t AliTRDCalibraFillHisto::UpdateHistograms(const AliTRDtrack *t) -{ - // - // Use AliTRDtrack for the calibration - // - - - AliTRDcluster *cl = 0x0; // pointeur to cluster - Int_t index0 = 0; // index of the first cluster in the current chamber - Int_t index1 = 0; // index of the current cluster in the current chamber - - ////////////////////////////// - // loop over the clusters - /////////////////////////////// - while((cl = t->GetCluster(index1))){ - - ///////////////////////////////////////////////////////////////////////// - // Fill the infos for the previous clusters if not the same detector - //////////////////////////////////////////////////////////////////////// - Int_t detector = cl->GetDetector(); - if ((detector != fDetectorPreviousTrack) && - (index0 != index1)) { - - fNumberTrack++; - - //If the same track, then look if the previous detector is in - //the same plane, if yes: not a good track - if ((GetLayer(detector) == GetLayer(fDetectorPreviousTrack))) { - return kFALSE; - } - - // Fill only if the track doesn't touch a masked pad or doesn't - if (fGoodTracklet) { - - // drift velocity unables to cut bad tracklets - Bool_t pass = FindP1TrackPHtracklet(t,index0,index1); - - // Gain calibration - if (fCH2dOn) { - FillTheInfoOfTheTrackCH(index1-index0); - } - - // PH calibration - if (fPH2dOn) { - FillTheInfoOfTheTrackPH(); - } - - if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1); - - - } // if a good tracklet - - // reset stuff - ResetfVariablestracklet(); - index0 = index1; - - } // Fill at the end the charge - - - ///////////////////////////////////////////////////////////// - // Localise and take the calib gain object for the detector - //////////////////////////////////////////////////////////// - if (detector != fDetectorPreviousTrack) { - - //Localise the detector - LocalisationDetectorXbins(detector); - - // Get cal - AliTRDcalibDB *cal = AliTRDcalibDB::Instance(); - if (!cal) { - AliInfo("Could not get calibDB"); - return kFALSE; - } - - // Get calib objects - if(!fIsHLT) InitCalPad(detector); - - } - - // Reset the detectbjobsor - fDetectorPreviousTrack = detector; - - /////////////////////////////////////// - // Store the info of the cluster - /////////////////////////////////////// - Int_t row = cl->GetPadRow(); - Int_t col = cl->GetPadCol(); - CheckGoodTrackletV1(cl); - Int_t group[2] = {0,0}; - if(fCH2dOn) group[0] = CalculateCalibrationGroup(0,row,col); - if(fPH2dOn) group[1] = CalculateCalibrationGroup(1,row,col); - StoreInfoCHPHtrack(cl,t->GetClusterdQdl(index1),group,row,col); - - index1++; - - } // while on clusters - - /////////////////////////// - // Fill the last plane - ////////////////////////// - if( index0 != index1 ){ - - fNumberTrack++; - - if (fGoodTracklet) { - - // drift velocity unables to cut bad tracklets - Bool_t pass = FindP1TrackPHtracklet(t,index0,index1); - - // Gain calibration - if (fCH2dOn) { - FillTheInfoOfTheTrackCH(index1-index0); - } - - // PH calibration - if (fPH2dOn) { - FillTheInfoOfTheTrackPH(); - } - - if(pass && fPRF2dOn) HandlePRFtracklet(t,index0,index1); - - } // if a good tracklet - - } - - // reset stuff - ResetfVariablestracklet(); - - return kTRUE; - -} -//____________Offline tracking in the AliTRDtracker____________________________ -Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) +Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t,const AliESDtrack *esdtrack) { // // Use AliTRDtrackV1 for the calibration @@ -754,15 +677,15 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) AliTRDcluster *cl = 0x0; // cluster attached now to the tracklet AliTRDcluster *cls = 0x0; // shared cluster attached now to the tracklet Bool_t newtr = kTRUE; // new track - - // Get cal - // AliTRDcalibDB *cal = AliTRDcalibDB::Instance(); - /* - if (!cal) { - AliInfo("Could not get calibDB"); - return kFALSE; - } -*/ + + + // + // Cut on the number of TRD tracklets + // + Int_t numberoftrdtracklets = t->GetNumberOfTracklets(); + if(numberoftrdtracklets < fMinNbTRDtracklets) return kFALSE; + + // if (!fCalibDB) { AliInfo("Could not get calibDB"); return kFALSE; @@ -778,6 +701,9 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) if(!tracklet->IsOK()) continue; fNumberTrack++; ResetfVariablestracklet(); + Float_t momentum = t->GetMomentum(itr); + if(TMath::Abs(momentum) < fMinTRDMomentum) continue; + ////////////////////////////////////////// // localisation of the tracklet and dqdl @@ -808,6 +734,7 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) //////////////////////////// // loop over the clusters //////////////////////////// + Double_t chargeQ = 0.0; Int_t nbclusters = 0; for(int jc=0; jcGetClusters(jc))) continue; @@ -823,7 +750,7 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) // Add the charge if shared cluster cls = tracklet->GetClusters(jc+AliTRDseedV1::kNtb); // - StoreInfoCHPHtrack(cl, tracklet->GetdQdl(jc),group,row,col,cls); + chargeQ += StoreInfoCHPHtrack(cl, tracklet->GetQperTB(jc),group,row,col,cls); //tracklet->GetdQdl(jc) } //////////////////////////////////////// @@ -838,16 +765,74 @@ Bool_t AliTRDCalibraFillHisto::UpdateHistogramsV1(const AliTRDtrackV1 *t) // Gain calibration if (fCH2dOn) { - FillTheInfoOfTheTrackCH(nbclusters); + if(fCutWithVdriftCalib) { + if(pass) FillTheInfoOfTheTrackCH(nbclusters); + } else { + FillTheInfoOfTheTrackCH(nbclusters); + } } // PH calibration if (fPH2dOn) { - FillTheInfoOfTheTrackPH(); + if(fCutWithVdriftCalib) { + if(pass) FillTheInfoOfTheTrackPH(); + } + else { + FillTheInfoOfTheTrackPH(); + } } if(pass && fPRF2dOn) HandlePRFtrackletV1(tracklet,nbclusters); - + + + ///////////////////////////////////////////////////////// + // Debug + //////////////////////////////////////////////////////// + if(fDebugLevel > 0){ + //printf("test\n"); + if ( !fDebugStreamer ) { + //debug stream + TDirectory *backup = gDirectory; + fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root"); + if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer + } + + Int_t stacke = AliTRDgeometry::GetStack(detector); + Int_t sme = AliTRDgeometry::GetSector(detector); + Int_t layere = AliTRDgeometry::GetLayer(detector); + // Some variables + Float_t b[2] = {0.0,0.0}; + Float_t bCov[3] = {0.0,0.0,0.0}; + if(esdtrack) esdtrack->GetImpactParameters(b,bCov); + if (bCov[0]<=0 || bCov[2]<=0) { + bCov[0]=0; bCov[2]=0; + } + Float_t dcaxy = b[0]; + Float_t dcaz = b[1]; + Int_t tpcnbclusters = 0; + if(esdtrack) tpcnbclusters = esdtrack->GetTPCclusters(0); + Double_t tpcsignal = 0.0; + if(esdtrack) tpcsignal = esdtrack->GetTPCsignal(); + Int_t cutvdriftlinear = 0; + if(!pass) cutvdriftlinear = 1; + + (* fDebugStreamer) << "FillCharge"<< + "detector="<GetCluster(index0))->GetDetector(); //detector - // parameters of the track - Double_t snp = t->GetSnpPlane(GetLayer(detector)); // sin angle in the plan yx track - Double_t tgl = t->GetTglPlane(GetLayer(detector)); // dz/dl and not dz/dx! - Double_t tnp = 0.0; // tan angle in the plan xy track - if( TMath::Abs(snp) < 1.){ - tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp)); - } - Float_t dzdx = tgl*TMath::Sqrt(1+tnp*tnp); // dz/dx now from dz/dl - // tilting pad and cross row - Int_t crossrow = 0; // if it crosses a pad row - Int_t rowp = -1; // if it crosses a pad row - AliTRDpadPlane *padplane = fGeo->GetPadPlane(GetLayer(detector),GetStack(detector)); - Double_t tiltingangle = padplane->GetTiltingAngle(); // tiltingangle of the pad - Float_t tnt = TMath::Tan(tiltingangle/180.*TMath::Pi()); // tan tiltingangle - // linear fit - fLinearFitterTracklet->ClearPoints(); - Double_t dydt = 0.0; // dydt tracklet after straight line fit - Double_t errorpar = 0.0; // error after straight line fit on dy/dt - Double_t pointError = 0.0; // error after straight line fit - Int_t nbli = 0; // number linear fitter points - - ////////////////////////////// - // loop over clusters - //////////////////////////// - for(Int_t k = 0; k < npoints; k++){ - - AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0); - if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue; - Double_t ycluster = cl->GetY(); - Int_t time = cl->GetPadTime(); - Double_t timeis = time/fSf; - //Double_t q = cl->GetQ(); - //See if cross two pad rows - Int_t row = cl->GetPadRow(); - if(k==0) rowp = row; - if(row != rowp) crossrow = 1; - - fLinearFitterTracklet->AddPoint(&timeis,ycluster,1); - nbli++; - - } - - ////////////////////////////// - // linear fit - ////////////////////////////// - if(nbli <= 2){ - fLinearFitterTracklet->ClearPoints(); - return kFALSE; - } - TVectorD pars; - fLinearFitterTracklet->Eval(); - fLinearFitterTracklet->GetParameters(pars); - pointError = TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2)); - errorpar = fLinearFitterTracklet->GetParError(1)*pointError; - dydt = pars[1]; - fLinearFitterTracklet->ClearPoints(); - - ///////////////////////////// - // debug - //////////////////////////// - if(fDebugLevel > 0){ - if ( !fDebugStreamer ) { - //debug stream - TDirectory *backup = gDirectory; - fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root"); - if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer - } - - - (* fDebugStreamer) << "FindP1TrackPHtracklet0"<< - //"snpright="< fNumberClustersf) return kFALSE; - if(pointError >= 0.1) return kFALSE; - if(crossrow == 1) return kFALSE; - - //////////////////////////// - // fill - //////////////////////////// - if(fLinearFitterOn){ - //Add to the linear fitter of the detector - if( TMath::Abs(snp) < 1.){ - Double_t x = tnp-dzdx*tnt; - if(fLinearFitterDebugOn) { - (GetLinearFitter(detector,kTRUE))->AddPoint(&x,dydt); - fEntriesLinearFitter[detector]++; - } - fLinearVdriftFit->Update(detector,x,pars[1]); - } - } - - return kTRUE; -} -//____________Offine tracking in the AliTRDtracker_____________________________ Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters) { // @@ -1056,414 +901,156 @@ Bool_t AliTRDCalibraFillHisto::FindP1TrackPHtrackletV1(const AliTRDseedV1 *track for(int icc=AliTRDseedV1::kNtb; iccGetClusters(icc); if(cl) crossrow = 1; - } - ////////////////////////////////// - // Loop clusters - ////////////////////////////////// - for(int ic=0; icGetClusters(ic))) continue; - if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue; - - Double_t ycluster = cl->GetY(); - Int_t time = cl->GetPadTime(); - Double_t timeis = time/fSf; - //See if cross two pad rows - Int_t row = cl->GetPadRow(); - if(rowp==-1) rowp = row; - if(row != rowp) crossrow = 1; - - fLinearFitterTracklet->AddPoint(&timeis,ycluster,1); - nbli++; - - - } - - //////////////////////////////////// - // Do the straight line fit now - /////////////////////////////////// - if(nbli <= 2){ - fLinearFitterTracklet->ClearPoints(); - return kFALSE; - } - TVectorD pars; - fLinearFitterTracklet->Eval(); - fLinearFitterTracklet->GetParameters(pars); - pointError = TMath::Sqrt(fLinearFitterTracklet->GetChisquare()/(nbli-2)); - errorpar = fLinearFitterTracklet->GetParError(1)*pointError; - dydt = pars[1]; - //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",fLinearFitterTracklet->GetChisquare(),nbli,pointError,fLinearFitterTracklet->GetParError(1),errorpar); - fLinearFitterTracklet->ClearPoints(); - - //////////////////////////////// - // Debug stuff - /////////////////////////////// - - - if(fDebugLevel > 0){ - if ( !fDebugStreamer ) { - //debug stream - TDirectory *backup = gDirectory; - fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root"); - if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer - } - - - Int_t layer = GetLayer(fDetectorPreviousTrack); - - (* fDebugStreamer) << "FindP1TrackPHtrackletV1"<< - //"snpright="< fNumberClustersf) return kFALSE; - if(pointError >= 0.3) return kFALSE; - if(crossrow == 1) return kFALSE; - - /////////////////////// - // Fill - ////////////////////// - - if(fLinearFitterOn){ - //Add to the linear fitter of the detector - if( TMath::Abs(snp) < 1.){ - Double_t x = tnp-dzdx*tnt; - if(fLinearFitterDebugOn) { - (GetLinearFitter(fDetectorPreviousTrack,kTRUE))->AddPoint(&x,dydt); - fEntriesLinearFitter[fDetectorPreviousTrack]++; - } - fLinearVdriftFit->Update(fDetectorPreviousTrack,x,pars[1]); - } - } - - return kTRUE; -} -//____________Offine tracking in the AliTRDtracker_____________________________ -Bool_t AliTRDCalibraFillHisto::HandlePRFtracklet(const AliTRDtrack *t, Int_t index0, Int_t index1) -{ - // - // PRF width calibration - // Assume a Gaussian shape: determinate the position of the three pad clusters - // Fit with a straight line - // Take the fitted values for all the clusters (3 or 2 pad clusters) - // Fill the PRF as function of angle of the track - // - // - - - ////////////////////////// - // variables - ///////////////////////// - Int_t npoints = index1-index0; // number of total points - Int_t nb3pc = 0; // number of three pads clusters used for fit - Int_t detector = ((AliTRDcluster *) t->GetCluster(index0))->GetDetector(); // detector - // To see the difference due to the fit - Double_t *padPositions; - padPositions = new Double_t[npoints]; - for(Int_t k = 0; k < npoints; k++){ - padPositions[k] = 0.0; - } - // Take the tgl and snp with the track t now - Double_t tgl = t->GetTglPlane(GetLayer(detector)); //dz/dl and not dz/dx - Double_t snp = t->GetSnpPlane(GetLayer(detector)); // sin angle in xy plan - Float_t dzdx = 0.0; // dzdx - Float_t tnp = 0.0; - if(TMath::Abs(snp) < 1.0){ - tnp = snp / TMath::Sqrt((1.-snp)*(1.+snp)); - dzdx = tgl*TMath::Sqrt(1+tnp*tnp); - } - // linear fitter - fLinearFitterTracklet->ClearPoints(); - - /////////////////////////// - // calcul the tnp group - /////////////////////////// - Bool_t echec = kFALSE; - Double_t shift = 0.0; - //Calculate the shift in x coresponding to this tnp - if(fNgroupprf != 0.0){ - shift = -3.0*(fNgroupprf-1)-1.5; - Double_t limithigh = -0.2*(fNgroupprf-1); - if((tnp < (-0.2*fNgroupprf)) || (tnp > (0.2*fNgroupprf))) echec = kTRUE; - else{ - while(tnp > limithigh){ - limithigh += 0.2; - shift += 3.0; - } - } - } - if(echec) { - delete [] padPositions; - return kFALSE; - } - - ////////////////////// - // loop clusters - ///////////////////// - for(Int_t k = 0; k < npoints; k++){ - //Take the cluster - AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0); - if(!cl) continue; - Short_t *signals = cl->GetSignals(); - Double_t time = cl->GetPadTime(); - //Calculate x if possible - Float_t xcenter = 0.0; - Bool_t echec1 = kTRUE; - if((time<=7) || (time>=21)) continue; - // Center 3 balanced: position with the center of the pad - if ((((Float_t) signals[3]) > 0.0) && - (((Float_t) signals[2]) > 0.0) && - (((Float_t) signals[4]) > 0.0)) { - echec1 = kFALSE; - // Security if the denomiateur is 0 - if ((((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) / - ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4])))) != 1.0) { - xcenter = 0.5 * (TMath::Log((Float_t) (((Float_t) signals[4]) / ((Float_t) signals[2])))) - / (TMath::Log(((Float_t) (((Float_t) signals[3]) * ((Float_t) signals[3]))) - / ((Float_t) (((Float_t) signals[2]) * ((Float_t) signals[4]))))); - } - else { - echec1 = kTRUE; - } - } - if(TMath::Abs(xcenter) > 0.5) echec = kTRUE; - if(echec) continue; - //if no echec: calculate with the position of the pad - // Position of the cluster - Double_t padPosition = xcenter + cl->GetPadCol(); - padPositions[k] = padPosition; - nb3pc++; - fLinearFitterTracklet->AddPoint(&time, padPosition,1); - }//clusters loop - - - ///////////////////////////// - // fit - //////////////////////////// - if(nb3pc < 3) { - delete [] padPositions; - fLinearFitterTracklet->ClearPoints(); - return kFALSE; - } - fLinearFitterTracklet->Eval(); - TVectorD line(2); - fLinearFitterTracklet->GetParameters(line); - Float_t pointError = -1.0; - if( fLinearFitterTracklet->GetChisquare()>=0.0) { - pointError = TMath::Sqrt( fLinearFitterTracklet->GetChisquare()/(nb3pc-2)); - } - fLinearFitterTracklet->ClearPoints(); - - ///////////////////////////////////////////////////// - // Now fill the PRF: second loop over clusters - ///////////////////////////////////////////////////// - for(Int_t k = 0; k < npoints; k++){ - //Take the cluster - AliTRDcluster *cl = (AliTRDcluster *) t->GetCluster(k+index0); - Short_t *signals = cl->GetSignals(); // signal - Double_t time = cl->GetPadTime(); // time bin - Float_t padPosTracklet = line[0]+line[1]*time; // reconstruct position from fit - Float_t padPos = cl->GetPadCol(); // middle pad - Double_t dpad = padPosTracklet - padPos; // reconstruct position relative to middle pad from fit - Float_t ycenter = 0.0; // relative center charge - Float_t ymin = 0.0; // relative left charge - Float_t ymax = 0.0; // relative right charge - - //Requiere simply two pads clusters at least - if(((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[2]) > 0.0)) || - ((((Float_t) signals[3]) > 0.0) && (((Float_t) signals[4]) > 0.0))){ - Float_t sum = ((Float_t) signals[2]) + ((Float_t) signals[3]) + ((Float_t) signals[4]); - if(sum > 0.0) ycenter = ((Float_t) signals[3])/ sum; - if(sum > 0.0) ymin = ((Float_t) signals[2])/ sum; - if(sum > 0.0) ymax = ((Float_t) signals[4])/ sum; - } - - //calibration group - Int_t rowcl = cl->GetPadRow(); // row of cluster - Int_t colcl = cl->GetPadCol(); // col of cluster - Int_t grouplocal = CalculateCalibrationGroup(2,rowcl,colcl); // calcul the corresponding group - Int_t caligroup = fCalibraMode->GetXbins(2)+ grouplocal; // calcul the corresponding group - Float_t xcl = cl->GetY(); // y cluster - Float_t qcl = cl->GetQ(); // charge cluster - Int_t layer = GetLayer(detector); // layer - Int_t stack = GetStack(detector); // stack - Double_t xdiff = dpad; // reconstructed position constant - Double_t x = dpad; // reconstructed position moved - Float_t ep = pointError; // error of fit - Float_t signal1 = (Float_t)signals[1]; // signal at the border - Float_t signal3 = (Float_t)signals[3]; // signal - Float_t signal2 = (Float_t)signals[2]; // signal - Float_t signal4 = (Float_t)signals[4]; // signal - Float_t signal5 = (Float_t)signals[5]; // signal at the border - - ////////////////////////////// - // debug - ///////////////////////////// - if(fDebugLevel > 0){ - if ( !fDebugStreamer ) { - //debug stream - TDirectory *backup = gDirectory; - fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root"); - if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer - } - - - x = xdiff; - Int_t type=0; - Float_t y = ycenter; - (* fDebugStreamer) << "HandlePRFtracklet"<< - "caligroup="< fNumberClustersf) continue; - if(nb3pc <= 5) continue; - if((time >= 21) || (time < 7)) continue; - if(TMath::Abs(snp) >= 1.0) continue; - if(TMath::Abs(qcl) < 80) continue; - - //////////////////////////// - // Fill - /////////////////////////// - if (fHisto2d) { - if(TMath::Abs(dpad) < 1.5) { - fPRF2d->Fill(shift+dpad,(caligroup+0.5),ycenter); - fPRF2d->Fill(shift-dpad,(caligroup+0.5),ycenter); - } - if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) { - fPRF2d->Fill(shift-(dpad+1.0),(caligroup+0.5),ymin); - fPRF2d->Fill(shift+(dpad+1.0),(caligroup+0.5),ymin); - } - if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) { - fPRF2d->Fill(shift+1.0-dpad,(caligroup+0.5),ymax); - fPRF2d->Fill(shift-1.0+dpad,(caligroup+0.5),ymax); - } + Float_t sigArr[AliTRDfeeParam::GetNcol()]; + memset(sigArr, 0, AliTRDfeeParam::GetNcol()*sizeof(sigArr[0])); + Int_t ncl=0, tbf=0, tbl=0; + + for(int ic=0; icGetClusters(ic))) continue; + + if(!tbf) tbf=ic; + tbl=ic; + ncl++; + Int_t col = cl->GetPadCol(); + for(int ip=-1, jp=2; jp<5; ip++, jp++){ + Int_t idx=col+ip; + if(idx<0 || idx>=AliTRDfeeParam::GetNcol()) continue; + sigArr[idx]+=((Float_t)cl->GetSignals()[jp]); } - if (fVector2d) { - if(TMath::Abs(dpad) < 1.5) { - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+dpad,ycenter); - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-dpad,ycenter); - } - if((ymin > 0.0) && (TMath::Abs(dpad+1.0) < 1.5)) { - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-(dpad+1.0),ymin); - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+(dpad+1.0),ymin); - } - if((ymax > 0.0) && (TMath::Abs(dpad-1.0) < 1.5)) { - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift+1.0-dpad,ymax); - fCalibraVector->UpdateVectorPRF(fDetectorPreviousTrack,grouplocal,shift-1.0+dpad,ymax); + + if((fLimitChargeIntegration) && (!cl->IsInChamber())) continue; + + Double_t ycluster = cl->GetY(); + Int_t time = cl->GetPadTime(); + Double_t timeis = time/fSf; + //See if cross two pad rows + Int_t row = cl->GetPadRow(); + if(rowp==-1) rowp = row; + if(row != rowp) crossrow = 1; + + fLinearFitterTracklet->AddPoint(&timeis,ycluster,1); + nbli++; + + + } + + //////////////////////////////////// + // Do the straight line fit now + /////////////////////////////////// + if(nbli <= 2){ + fLinearFitterTracklet->ClearPoints(); + return kFALSE; + } + TVectorD pars; + fLinearFitterTracklet->Eval(); + fLinearFitterTracklet->GetParameters(pars); + pointError = TMath::Sqrt(TMath::Abs(fLinearFitterTracklet->GetChisquare()/(nbli-2))); + errorpar = fLinearFitterTracklet->GetParError(1)*pointError; + dydt = pars[1]; + //printf("chis %f, nbli %d, pointError %f, parError %f, errorpar %f\n",fLinearFitterTracklet->GetChisquare(),nbli,pointError,fLinearFitterTracklet->GetParError(1),errorpar); + fLinearFitterTracklet->ClearPoints(); + + //////////////////////////////////// + // Calc the projection of the clusters on the y direction + /////////////////////////////////// + + Float_t signalSum(0.); + Float_t mean = 0.0, rms = 0.0; + Float_t dydx = tracklet->GetYref(1), tilt = tracklet->GetTilt(); // ,dzdx = tracklet->GetZref(1); (identical to the previous definition!) + Float_t dz = dzdx*(tbl-tbf)/10; + if(ncl>10){ + for(Int_t ip(0); ip 0.0) mean/=signalSum; + + for(Int_t ip = 0; ip 0.0) rms = TMath::Sqrt(TMath::Abs(rms/signalSum)); + + rms -= TMath::Abs(dz*tilt); + dydx -= dzdx*tilt; + } + + //////////////////////////////// + // Debug stuff + /////////////////////////////// + + + if(fDebugLevel > 1){ + if ( !fDebugStreamer ) { + //debug stream + TDirectory *backup = gDirectory; + fDebugStreamer = new TTreeSRedirector("TRDdebugCalibraFill.root"); + if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer + } + + float xcoord = tnp-dzdx*tnt; + float pt = tracklet->GetPt(); + Int_t layer = GetLayer(fDetectorPreviousTrack); + + (* fDebugStreamer) << "FindP1TrackPHtrackletV1"<< + //"snpright="< fNumberClustersf) return kFALSE; + if(pointError >= 0.3) return kFALSE; + if(crossrow == 1) return kTRUE; + + /////////////////////// + // Fill + ////////////////////// + + if(fLinearFitterOn){ + //Add to the linear fitter of the detector + if( TMath::Abs(snp) < 1.){ + Double_t x = tnp-dzdx*tnt; + if(fLinearFitterDebugOn) { + (GetLinearFitter(fDetectorPreviousTrack,kTRUE))->AddPoint(&x,dydt); + fEntriesLinearFitter[fDetectorPreviousTrack]++; } + fLinearVdriftFit->Update(fDetectorPreviousTrack,x,pars[1]); } } - delete [] padPositions; - return kTRUE; + if(fExbAltFitOn){ + fExbAltFit->Update(fDetectorPreviousTrack,dydx,rms); + } + return kTRUE; } //____________Offine tracking in the AliTRDtracker_____________________________ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet, Int_t nbclusters) @@ -1656,7 +1243,7 @@ Bool_t AliTRDCalibraFillHisto::HandlePRFtrackletV1(const AliTRDseedV1 *tracklet, // Debug stuff //////////////////// - if(fDebugLevel > 0){ + if(fDebugLevel > 1){ if ( !fDebugStreamer ) { //debug stream TDirectory *backup = gDirectory; @@ -1859,8 +1446,8 @@ Bool_t AliTRDCalibraFillHisto::IsPadOn(Int_t detector, Int_t row, Int_t col) con return kFALSE; } - if (!cal->IsChamberInstalled(detector) || - cal->IsChamberMasked(detector) || + if (!cal->IsChamberGood(detector) || + cal->IsChamberNoData(detector) || cal->IsPadMasked(detector,col,row)) { return kFALSE; } @@ -2042,12 +1629,14 @@ void AliTRDCalibraFillHisto::SetNumberGroupsPRF(Short_t numberGroupsPRF) // Per tracklet: store or reset the info, fill the histos with the info ////////////////////////////////////////////////////////////////////////////////////////// //_____________________________________________________________________________ -void AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Double_t dqdl,const Int_t *group,const Int_t row,const Int_t col,const AliTRDcluster *cls) +Float_t AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Double_t dqdl,const Int_t *group,const Int_t row,const Int_t col,const AliTRDcluster *cls) { // // Store the infos in fAmpTotal, fPHPlace and fPHValue // Correct from the gain correction before // cls is shared cluster if any + // Return the charge + // // //printf("StoreInfoCHPHtrack\n"); @@ -2067,7 +1656,7 @@ void AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Do if(fIsHLT) correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack); else correctthegain = fCalDetGain->GetValue(fDetectorPreviousTrack)*fCalROCGain->GetValue(col,row); Float_t correction = 1.0; - Float_t normalisation = 6.67; + Float_t normalisation = 1.13; //org: 6.67; 1st: 1.056; 2nd: 1.13; // we divide with gain in AliTRDclusterizer::Transform... if( correctthegain > 0 ) normalisation /= correctthegain; @@ -2089,6 +1678,8 @@ void AliTRDCalibraFillHisto::StoreInfoCHPHtrack(const AliTRDcluster *cl,const Do fPHPlace[time] = group[1]; fPHValue[time] = charge; } + + return correction; } //____________Offine tracking in the AliTRDtracker_____________________________ @@ -2229,7 +1820,7 @@ void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackCH(Int_t nbclusters) } break; default: break; - } + } } //____________Offine tracking in the AliTRDtracker_____________________________ void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackPH() @@ -2442,227 +2033,7 @@ void AliTRDCalibraFillHisto::FillTheInfoOfTheTrackPH() // DAQ process functions ///////////////////////////////////////////////////////////////////////////////////////// //_____________________________________________________________________ -Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliTRDrawStreamBase *rawStream, Bool_t nocheck) -{ - // - // Event Processing loop - AliTRDrawStreamBase - // TestBeam 2007 version - // 0 timebin problem - // 1 no input - // 2 input - // Same algorithm as TestBeam but different reader - // - - rawStream->SetSharedPadReadout(kFALSE); - - Int_t withInput = 1; - - Double_t phvalue[16][144][36]; - for(Int_t k = 0; k < 36; k++){ - for(Int_t j = 0; j < 16; j++){ - for(Int_t c = 0; c < 144; c++){ - phvalue[j][c][k] = 0.0; - } - } - } - - fDetectorPreviousTrack = -1; - fMCMPrevious = -1; - fROBPrevious = -1; - Int_t nbtimebin = 0; - Int_t baseline = 10; - //printf("------------Detector\n"); - - if(!nocheck){ - - fTimeMax = 0; - - while (rawStream->Next()) { - - Int_t idetector = rawStream->GetDet(); // current detector - Int_t imcm = rawStream->GetMCM(); // current MCM - Int_t irob = rawStream->GetROB(); // current ROB - - //printf("Detector %d\n",idetector); - - if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){ - - // Fill - withInput = TMath::Max(FillDAQ(phvalue),withInput); - - - // reset - for(Int_t k = 0; k < 36; k++){ - for(Int_t j = 0; j < 16; j++){ - for(Int_t c = 0; c < 144; c++){ - phvalue[j][c][k] = 0.0; - } - } - } - } - - fDetectorPreviousTrack = idetector; - fMCMPrevious = imcm; - fROBPrevious = irob; - - nbtimebin = rawStream->GetNumberOfTimeBins(); // number of time bins read from data - if(nbtimebin == 0) return 0; - if((fTimeMax != 0) && (nbtimebin != fTimeMax)) return 0; - fTimeMax = nbtimebin; - - //baseline = rawStream->GetCommonAdditive(); // common additive baseline - fNumberClustersf = fTimeMax; - fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax); - - - Int_t *signal = rawStream->GetSignals(); // current ADC signal - Int_t col = rawStream->GetCol(); - Int_t row = rawStream->GetRow(); - - - // printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline); - - - for(Int_t itime = 0; itime < nbtimebin; itime++){ - phvalue[row][col][itime] = signal[itime]-baseline; - } - } - - // fill the last one - if(fDetectorPreviousTrack != -1){ - - // Fill - withInput = TMath::Max(FillDAQ(phvalue),withInput); - - // reset - for(Int_t k = 0; k < 36; k++){ - for(Int_t j = 0; j < 16; j++){ - for(Int_t c = 0; c < 144; c++){ - phvalue[j][c][k] = 0.0; - } - } - } - } - - } - else{ - - while (rawStream->Next()) { //iddetecte - - Int_t idetector = rawStream->GetDet(); // current detector - Int_t imcm = rawStream->GetMCM(); // current MCM - Int_t irob = rawStream->GetROB(); // current ROB - - //printf("Detector %d\n",idetector); - - if((fDetectorPreviousTrack != idetector) && (fDetectorPreviousTrack != -1)){ - - // Fill - withInput = TMath::Max(FillDAQ(phvalue),withInput); - - // reset - for(Int_t k = 0; k < 36; k++){ - for(Int_t j = 0; j < 16; j++){ - for(Int_t c = 0; c < 144; c++){ - phvalue[j][c][k] = 0.0; - } - } - } - } - - fDetectorPreviousTrack = idetector; - fMCMPrevious = imcm; - fROBPrevious = irob; - - //baseline = rawStream->GetCommonAdditive(); // common baseline - - fTimeMax = rawStream->GetNumberOfTimeBins(); // number of time bins read from data - fNumberClustersf = fTimeMax; - fNumberClusters = (Int_t)(fNumberClustersProcent*fTimeMax); - Int_t *signal = rawStream->GetSignals(); // current ADC signal - Int_t col = rawStream->GetCol(); - Int_t row = rawStream->GetRow(); - - - //printf("detector %d, signal[0] %d, signal[1] %d, signal[2] %d, baseline %d\n",idetector,signal[0],signal[1],signal[2], baseline); - - for(Int_t itime = 0; itime < fTimeMax; itime++){ - phvalue[row][col][itime] = signal[itime]-baseline; - /*if(phvalue[row][col][itime] >= 20) { - printf("----------> phvalue[%d][%d][%d] %d baseline %d \n", - row, - col, - itime, - signal[itime], - baseline); - }*/ - } - } - - // fill the last one - if(fDetectorPreviousTrack != -1){ - - // Fill - withInput = TMath::Max(FillDAQ(phvalue),withInput); - - // reset - for(Int_t k = 0; k < 36; k++){ - for(Int_t j = 0; j < 16; j++){ - for(Int_t c = 0; c < 144; c++){ - phvalue[j][c][k] = 0.0; - } - } - } - } - } - - return withInput; - -} -//_____________________________________________________________________ -Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliRawReader *rawReader, Bool_t nocheck) -{ - // - // Event processing loop - AliRawReader - // Testbeam 2007 version - // - - AliTRDrawStreamBase rawStream(rawReader); - - rawReader->Select("TRD"); - - return ProcessEventDAQ(&rawStream, nocheck); -} - -//_________________________________________________________________________ -Int_t AliTRDCalibraFillHisto::ProcessEventDAQ( -#ifdef ALI_DATE - const eventHeaderStruct *event, - Bool_t nocheck -#else - const eventHeaderStruct* /*event*/, - Bool_t /*nocheck*/ - -#endif - ) -{ - // - // process date event - // Testbeam 2007 version - // -#ifdef ALI_DATE - AliRawReader *rawReader = new AliRawReaderDate((void*)event); - Int_t result=ProcessEventDAQ(rawReader, nocheck); - delete rawReader; - return result; -#else - Fatal("AliTRDCalibraFillHisto", "this class was compiled without DATE"); - return 0; -#endif - -} -//_____________________________________________________________________ -Int_t AliTRDCalibraFillHisto::ProcessEventDAQ2(AliRawReader *rawReader) +Int_t AliTRDCalibraFillHisto::ProcessEventDAQ(AliRawReader *rawReader) { //main // // Event Processing loop - AliTRDrawStream @@ -2674,8 +2045,6 @@ Int_t AliTRDCalibraFillHisto::ProcessEventDAQ2(AliRawReader *rawReader) // AliTRDrawStream *rawStream = new AliTRDrawStream(rawReader); - rawStream->SetNoErrorWarning(); - rawStream->SetSharedPadReadout(kFALSE); AliTRDdigitsManager *digitsManager = new AliTRDdigitsManager(kTRUE); digitsManager->CreateArrays(); @@ -2906,7 +2275,7 @@ Int_t AliTRDCalibraFillHisto::FillDAQ(Double_t phvalue[16][144][36]){ ///////////////////////////////////////////////////////// // Debug //////////////////////////////////////////////////////// - if(fDebugLevel > 0){ + if(fDebugLevel > 1){ if ( !fDebugStreamer ) { //debug stream TDirectory *backup = gDirectory; @@ -3310,6 +2679,8 @@ void AliTRDCalibraFillHisto::CreatePH2d(Int_t nn) name += fVersionVdriftUsed; name += "Subver"; name += fSubVersionVdriftUsed; + name += "FirstRun"; + name += fFirstRunVdrift; name += "Nz"; name += fCalibraMode->GetNz(1); name += "Nrphi"; @@ -3335,6 +2706,8 @@ void AliTRDCalibraFillHisto::CreateCH2d(Int_t nn) name += fVersionGainUsed; name += "Subver"; name += fSubVersionGainUsed; + name += "FirstRun"; + name += fFirstRunGain; name += "Nz"; name += fCalibraMode->GetNz(0); name += "Nrphi"; @@ -3499,15 +2872,20 @@ void AliTRDCalibraFillHisto::AnalyseLinearFitter() TVectorD *par = new TVectorD(2); TVectorD pare = TVectorD(2); TVectorD *parE = new TVectorD(3); - linearfitter->Eval(); - linearfitter->GetParameters(*par); - linearfitter->GetErrors(pare); - Float_t ppointError = TMath::Sqrt(TMath::Abs(linearfitter->GetChisquare())/fEntriesLinearFitter[k]); - (*parE)[0] = pare[0]*ppointError; - (*parE)[1] = pare[1]*ppointError; - (*parE)[2] = (Double_t) fEntriesLinearFitter[k]; - ((TObjArray *)fLinearVdriftFit->GetPArray())->AddAt(par,k); - ((TObjArray *)fLinearVdriftFit->GetEArray())->AddAt(parE,k); + if((linearfitter->EvalRobust(0.8)==0)) { + //linearfitter->Eval(); + linearfitter->GetParameters(*par); + //linearfitter->GetErrors(pare); + //Float_t ppointError = TMath::Sqrt(TMath::Abs(linearfitter->GetChisquare())/fEntriesLinearFitter[k]); + //(*parE)[0] = pare[0]*ppointError; + //(*parE)[1] = pare[1]*ppointError; + + (*parE)[0] = 0.0; + (*parE)[1] = 0.0; + (*parE)[2] = (Double_t) fEntriesLinearFitter[k]; + ((TObjArray *)fLinearVdriftFit->GetPArray())->AddAt(par,k); + ((TObjArray *)fLinearVdriftFit->GetEArray())->AddAt(parE,k); + } } } }