#include "DSelector_kpi0sigma0.h" void DSelector_kpi0sigma0::Init(TTree *locTree) { // The Init() function is called when the selector needs to initialize a new tree or chain. // Typically here the branch addresses and branch pointers of the tree will be set. // Init() will be called many times when running on PROOF (once per file to be processed). //SET OUTPUT FILE NAME //can be overriden by user in PROOF dOutputFileName = "kpi0sigma0.root"; //"" for none dOutputTreeFileName = ""; //"" for none //DO THIS NEXT //Because this function gets called for each TTree in the TChain, we must be careful: //We need to re-initialize the tree interface & branch wrappers, but don't want to recreate histograms bool locInitializedPriorFlag = dInitializedFlag; //save whether have been initialized previously DSelector::Init(locTree); //This must be called to initialize wrappers for each new TTree //gDirectory now points to the output file with name dOutputFileName (if any) if(locInitializedPriorFlag) return; //have already created histograms, etc. below: exit //THEN THIS Get_ComboWrappers(); dPreviousRunNumber = 0; /*********************************** EXAMPLE USER INITIALIZATION: ANALYSIS ACTIONS **********************************/ //ANALYSIS ACTIONS: //Executed in order if added to dAnalysisActions //false/true below: use measured/kinfit data //PID dAnalysisActions.push_back(new DHistogramAction_ParticleID(dComboWrapper, false)); dAnalysisActions.push_back(new DCutAction_dEdx(dComboWrapper, false, KPlus, SYS_CDC)); dAnalysisActions.push_back(new DCutAction_dEdx(dComboWrapper, false, Proton, SYS_CDC)); dAnalysisActions.push_back(new DCutAction_dEdx(dComboWrapper, false, PiMinus, SYS_CDC)); //below: value: +/- N ns, Unknown: All PIDs, SYS_NULL: all timing systems dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.4, KPlus, SYS_BCAL)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.25, KPlus, SYS_TOF)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.5, KPlus, SYS_FCAL)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 1.0, Proton, SYS_FCAL)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.8, Proton, SYS_BCAL)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.5, Proton, SYS_TOF)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 0.4, Gamma, SYS_BCAL)); dAnalysisActions.push_back(new DCutAction_PIDDeltaT(dComboWrapper, false, 1.0, Gamma, SYS_FCAL)); dAnalysisActions.push_back(new DHistogramAction_ParticleID(dComboWrapper, false, "PostPID")); //MASSES deque locKStarPIDs; locKStarPIDs.push_back(KPlus); locKStarPIDs.push_back(Pi0); deque locYStar0PIDs; locYStar0PIDs.push_back(Pi0); locYStar0PIDs.push_back(Sigma0); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Pi0, 1020, 0.05, 0.22, "Pi0")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Lambda, 1000, 1.0, 1.2, "Lambda")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Sigma0, 1000, 1.1, 1.3, "Sigma0")); dAnalysisActions.push_back(new DHistogramAction_MissingMassSquared(dComboWrapper, false, 1000, -0.1, 0.1)); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0_KinFit")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus_KinFit")); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, false, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0)); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, true, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0, "KinFit")); //KINFIT RESULTS dAnalysisActions.push_back(new DHistogramAction_KinFitResults(dComboWrapper)); dAnalysisActions.push_back(new DCutAction_KinFitFOM(dComboWrapper, 1.0E-3)); //MASSES dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Pi0, 1020, 0.05, 0.22, "Pi0_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Lambda, 1000, 1.0, 1.2, "Lambda_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Sigma0, 1000, 1.1, 1.3, "Sigma0_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_MissingMassSquared(dComboWrapper, false, 1000, -0.1, 0.1, "PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0_KinFit_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus_KinFit_PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, false, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0, "PostKinFitCut")); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, true, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0, "KinFit_PostKinFitCut")); //BEAM ENERGY dAnalysisActions.push_back(new DHistogramAction_BeamEnergy(dComboWrapper, false)); dAnalysisActions.push_back(new DCutAction_BeamEnergy(dComboWrapper, false, 8.4, 9.05)); //MASSES dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Pi0, 1020, 0.05, 0.22, "Pi0_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Lambda, 1000, 1.0, 1.2, "Lambda_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, Sigma0, 1000, 1.1, 1.3, "Sigma0_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_MissingMassSquared(dComboWrapper, false, 1000, -0.1, 0.1, "PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locKStarPIDs, 600, 0.6, 1.2, "KStar0_KinFit_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, false, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_InvariantMass(dComboWrapper, true, 0, locYStar0PIDs, 1200, 1.2, 2.4, "YStarPlus_KinFit_PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, false, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0, "PostBeamECut")); dAnalysisActions.push_back(new DHistogramAction_Dalitz(dComboWrapper, true, 0, locKStarPIDs, locYStar0PIDs, 700, 0.0, 14.0, 850, 1.0, 18.0, "KinFit_PostBeamECut")); //KINEMATICS dAnalysisActions.push_back(new DHistogramAction_ParticleComboKinematics(dComboWrapper, false)); dAnalysisActions.push_back(new DHistogramAction_ParticleComboKinematics(dComboWrapper, true, "KinFit")); //INITIALIZE ACTIONS //If you create any actions that you want to run manually (i.e. don't add to dAnalysisActions), be sure to initialize them here as well Initialize_Actions(); /******************************** EXAMPLE USER INITIALIZATION: STAND-ALONE HISTOGRAMS *******************************/ //EXAMPLE MANUAL HISTOGRAMS: dHist_NumSurvivingBeamParticles = new TH1I("NumSurvivingBeamParticles", "# Surviving Beam Particles", 100, 0.5, 100.5); //BEST COMBO HISTOGRAMS gDirectory->mkdir("BestCombo", "BestCombo")->cd(); { dHist_KStarMass_BestCombo = new TH1I("KStarMass", ";K^{#plus}#pi^{#minus} Invariant Mass (GeV/c^{2})", 600, 0.6, 1.2); dHist_MandelstamT = new TH1I("MandelstamT", ";t (GeV^{2})", 200, -5.0, 0.0); dHist_KStarMassVsT = new TH2I("KStarMassVsT", ";t (GeV^{2});K^{#plus}#pi^{#minus} Invariant Mass (GeV/c^{2})", 200, -5.0, 0.0, 300, 0.6, 1.2); //BEST COMBO: POLARIZATION dHist_ProdPlanePhi_PARA = new TH1I("ProdPlanePhi_PARA", ";Production Plane #phi#circ", 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Sideband_PARA = new TH2I("ProdPlanePhiVsT_Sideband_PARA", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Signal_PARA = new TH2I("ProdPlanePhiVsT_Signal_PARA", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); dHist_ProdPlanePhi_PERP = new TH1I("ProdPlanePhi_PERP", ";Production Plane #phi#circ", 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Sideband_PERP = new TH2I("ProdPlanePhiVsT_Sideband_PERP", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Signal_PERP = new TH2I("ProdPlanePhiVsT_Signal_PERP", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); dHist_ProdPlanePhi_AMO = new TH1I("ProdPlanePhi_AMO", ";Production Plane #phi#circ", 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Sideband_AMO = new TH2I("ProdPlanePhiVsT_Sideband_AMO", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); dHist_ProdPlanePhiVsT_Signal_AMO = new TH2I("ProdPlanePhiVsT_Signal_AMO", ";t (GeV^{2});Production Plane #phi#circ", 200, -5.0, 0.0, 360, -180.0, 180.0); } gDirectory->cd(".."); /***************************************** ADVANCED: CHOOSE BRANCHES TO READ ****************************************/ //TO SAVE PROCESSING TIME //If you know you don't need all of the branches/data, but just a subset of it, you can speed things up //By default, for each event, the data is retrieved for all branches //If you know you only need data for some branches, you can skip grabbing data from the branches you don't need //Do this by doing something similar to the commented code below //dTreeInterface->Clear_GetEntryBranches(); //now get none //dTreeInterface->Register_GetEntryBranch("Proton__P4"); //manually set the branches you want } Bool_t DSelector_kpi0sigma0::Process(Long64_t locEntry) { // The Process() function is called for each entry in the tree. The entry argument // specifies which entry in the currently loaded tree is to be processed. // // This function should contain the "body" of the analysis. It can contain // simple or elaborate selection criteria, run algorithms on the data // of the event and typically fill histograms. // // The processing can be stopped by calling Abort(). // Use fStatus to set the return value of TTree::Process(). // The return value is currently not used. //CALL THIS FIRST DSelector::Process(locEntry); //Gets the data from the tree for the entry //cout << "RUN " << Get_RunNumber() << ", EVENT " << Get_EventNumber() << endl; /******************************************** GET POLARIZATION ORIENTATION ******************************************/ //Only if the run number changes //RCDB environment must be setup in order for this to work! (Will return false otherwise) UInt_t locRunNumber = Get_RunNumber(); if(locRunNumber != dPreviousRunNumber) { dIsPolarizedFlag = dAnalysisUtilities.Get_IsPolarizedBeam(locRunNumber, dIsPARAFlag); dPreviousRunNumber = locRunNumber; } /********************************************* SETUP UNIQUENESS TRACKING ********************************************/ //ANALYSIS ACTIONS: Reset uniqueness tracking for each action //For any actions that you are executing manually, be sure to call Reset_NewEvent() on them here Reset_Actions_NewEvent(); //PREVENT-DOUBLE COUNTING WHEN HISTOGRAMMING //Sometimes, some content is the exact same between one combo and the next //e.g. maybe two combos have different beam particles, but the same data for the final-state //When histogramming, you don't want to double-count when this happens: artificially inflates your signal (or background) //So, for each quantity you histogram, keep track of what particles you used (for a given combo) //Then for each combo, just compare to what you used before, and make sure it's unique //EXAMPLE 1: Particle-specific info: set locUsedSoFar_BeamEnergy; //Int_t: Unique ID for beam particles. set: easy to use, fast to search //EXAMPLE 2: Combo-specific info: //In general: Could have multiple particles with the same PID: Use a set of Int_t's //In general: Multiple PIDs, so multiple sets: Contain within a map //Multiple combos: Contain maps within a set (easier, faster to search) set > > locUsedSoFar_MissingMass; //INSERT USER ANALYSIS UNIQUENESS TRACKING HERE /************************************************* LOOP OVER COMBOS *************************************************/ //Loop over combos UInt_t locNumSurvivingCombos = 0; set locSurvivingBeamIDs; double locBestMissingMassSquared = 9.9E9; double locT_BestCombo = 0.0, locKStarMass_BestCombo = 0.0, locProdPlanePhi_BestCombo = 0.0; bool locIsKStarSignal_BestCombo = false; for(UInt_t loc_i = 0; loc_i < Get_NumCombos(); ++loc_i) { //Set branch array indices for combo and all combo particles dComboWrapper->Set_ComboIndex(loc_i); // Is used to indicate when combos have been cut if(dComboWrapper->Get_IsComboCut()) // Is false when tree originally created continue; // Combo has been cut previously /********************************************** GET PARTICLE INDICES *********************************************/ //Used for tracking uniqueness when filling histograms, and for determining unused particles //Step 0 Int_t locBeamID = dComboBeamWrapper->Get_BeamID(); Int_t locKPlusTrackID = dKPlusWrapper->Get_TrackID(); //Step 1 Int_t locPhoton1NeutralID = dPhoton1Wrapper->Get_NeutralID(); Int_t locPhoton2NeutralID = dPhoton2Wrapper->Get_NeutralID(); //Step 2 Int_t locPhoton3NeutralID = dPhoton3Wrapper->Get_NeutralID(); //Step 3 Int_t locProtonTrackID = dProtonWrapper->Get_TrackID(); Int_t locPiMinusTrackID = dPiMinusWrapper->Get_TrackID(); /*********************************************** GET FOUR-MOMENTUM **********************************************/ // Get P4's: //is kinfit if kinfit performed, else is measured //dTargetP4 is target p4 //Step 0 TLorentzVector locBeamP4 = dComboBeamWrapper->Get_P4(); TLorentzVector locKPlusP4 = dKPlusWrapper->Get_P4(); //Step 1 TLorentzVector locDecayingPi0P4 = dDecayingPi0Wrapper->Get_P4(); TLorentzVector locPhoton1P4 = dPhoton1Wrapper->Get_P4(); TLorentzVector locPhoton2P4 = dPhoton2Wrapper->Get_P4(); //Step 2 TLorentzVector locDecayingSigma0P4 = dDecayingSigma0Wrapper->Get_P4(); TLorentzVector locPhoton3P4 = dPhoton3Wrapper->Get_P4(); //Step 3 TLorentzVector locDecayingLambdaP4 = dDecayingLambdaWrapper->Get_P4(); TLorentzVector locProtonP4 = dProtonWrapper->Get_P4(); TLorentzVector locPiMinusP4 = dPiMinusWrapper->Get_P4(); // Get Measured P4's: //Step 0 TLorentzVector locBeamP4_Measured = dComboBeamWrapper->Get_P4_Measured(); TLorentzVector locKPlusP4_Measured = dKPlusWrapper->Get_P4_Measured(); //Step 1 TLorentzVector locPhoton1P4_Measured = dPhoton1Wrapper->Get_P4_Measured(); TLorentzVector locPhoton2P4_Measured = dPhoton2Wrapper->Get_P4_Measured(); //Step 2 TLorentzVector locPhoton3P4_Measured = dPhoton3Wrapper->Get_P4_Measured(); //Step 3 TLorentzVector locProtonP4_Measured = dProtonWrapper->Get_P4_Measured(); TLorentzVector locPiMinusP4_Measured = dPiMinusWrapper->Get_P4_Measured(); /********************************************* COMBINE FOUR-MOMENTUM ********************************************/ // DO YOUR STUFF HERE // Combine 4-vectors TLorentzVector locMissingP4_Measured = locBeamP4_Measured + dTargetP4; locMissingP4_Measured -= locKPlusP4_Measured + locPhoton1P4_Measured + locPhoton2P4_Measured + locPhoton3P4_Measured + locProtonP4_Measured + locPiMinusP4_Measured; TLorentzVector locKStarP4_Measured = locKPlusP4_Measured + locPhoton1P4_Measured + locPhoton2P4_Measured; TLorentzVector locKStarP4_KinFit = locKPlusP4 + locPhoton1P4 + locPhoton2P4; if(dPhoton3Wrapper->Get_Energy_FCAL() > 0.0) continue; //photon from sigma0 /******************************************** EXECUTE ANALYSIS ACTIONS *******************************************/ // Loop through the analysis actions, executing them in order for the active particle combo if(!Execute_Actions()) //if the active combo fails a cut, IsComboCutFlag automatically set continue; //if you manually execute any actions, and it fails a cut, be sure to call: //dComboWrapper->Set_IsComboCut(true); /************************************************* MISCELLANEOUS *************************************************/ ++locNumSurvivingCombos; locSurvivingBeamIDs.insert(locBeamID); double locKStarMass = locKStarP4_Measured.M(); double locKStarMass_KinFit = locKStarP4_KinFit.M(); bool locIsKStarSignalFlag = ((locKStarMass_KinFit > 0.8) && (locKStarMass_KinFit <= 0.98)); double locProdPlanePhi = 0.0; //dAnalysisUtilities.Calc_ProdPlanePhi_Pseudoscalar(locBeamP4.E(), dTargetPID, locKPlusP4); double locT = (locBeamP4 - locKStarP4_KinFit).Mag2(); /**************************************************** BEST COMBO ****************************************************/ // Signal selection and background subtraction are performed on the KStar mass peak. // However, a given KStar can be matched to multiple beam photons (or kaon candidates) // This is OK in a 1D plot (check for duplicates), but we are quoting results in 2D: In bins of Mandelstam-t // We cannot put a signal KStar mass value into two separate t-bins: We would count that signal event in both t-bins. // Therefore, we must choose a t-bin. This is, in effect, a choice of "Best" combo. // We must choose this "best" combo in a way that doesn't bias the background subtraction. // Since the background subtraction is done on the KStar mass peak, we must choose something else. // Ideally we would pick the combo with the best kinematic-fit confidence level (or BDT response when it's ready). // However, at the time of this writing, the tracking errors are bad and this quantity cannot be trusted. // Therefore, we choose the combo that has a missing mass squared closest to zero: It's calculation involves all particles. // As far as the code is concerned, it is easiest to keep track of the relevant quantities for the "Best" combo // and update them as new "Best" combos are found. Then histogram the "Best" results at the end. double locMissingMassSquared = locMissingP4_Measured.M2(); if(fabs(locMissingMassSquared) > fabs(locBestMissingMassSquared)) continue; locBestMissingMassSquared = locMissingMassSquared; locT_BestCombo = locT; locKStarMass_BestCombo = locKStarMass_KinFit; locProdPlanePhi_BestCombo = locProdPlanePhi; locIsKStarSignal_BestCombo = locIsKStarSignalFlag; } // end of combo loop //FILL HISTOGRAMS: Num combos / events surviving actions Fill_NumCombosSurvivedHists(); dHist_NumSurvivingBeamParticles->Fill(locSurvivingBeamIDs.size()); //FILL ASYMMETRY INFO FOR BEST COMBOS if(locNumSurvivingCombos > 0) { dHist_KStarMass_BestCombo->Fill(locKStarMass_BestCombo); dHist_KStarMassVsT->Fill(locT_BestCombo, locKStarMass_BestCombo); dHist_MandelstamT->Fill(locT_BestCombo); if(locIsKStarSignal_BestCombo) { if(!dIsPolarizedFlag) { dHist_ProdPlanePhiVsT_Signal_AMO->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); dHist_ProdPlanePhi_AMO->Fill(locProdPlanePhi_BestCombo); } else if(dIsPARAFlag) { dHist_ProdPlanePhiVsT_Signal_PARA->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); dHist_ProdPlanePhi_PARA->Fill(locProdPlanePhi_BestCombo); } else //PERP { dHist_ProdPlanePhiVsT_Signal_PERP->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); dHist_ProdPlanePhi_PERP->Fill(locProdPlanePhi_BestCombo); } } else //sideband { if(!dIsPolarizedFlag) dHist_ProdPlanePhiVsT_Sideband_AMO->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); else if(dIsPARAFlag) dHist_ProdPlanePhiVsT_Sideband_PARA->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); else //PERP dHist_ProdPlanePhiVsT_Sideband_PERP->Fill(locT_BestCombo, locProdPlanePhi_BestCombo); } } /************************************ EXAMPLE: FILL CLONE OF TTREE HERE WITH CUTS APPLIED ************************************/ /* Bool_t locIsEventCut = true; for(UInt_t loc_i = 0; loc_i < Get_NumCombos(); ++loc_i) { //Set branch array indices for combo and all combo particles dComboWrapper->Set_ComboIndex(loc_i); // Is used to indicate when combos have been cut if(dComboWrapper->Get_IsComboCut()) continue; locIsEventCut = false; // At least one combo succeeded break; } if(!locIsEventCut && dOutputTreeFileName != "") FillOutputTree(); */ return kTRUE; } void DSelector_kpi0sigma0::Finalize(void) { //Save anything to output here that you do not want to be in the default DSelector output ROOT file. //Otherwise, don't do anything else (especially if you are using PROOF). //If you are using PROOF, this function is called on each thread, //so anything you do will not have the combined information from the various threads. //Besides, it is best-practice to do post-processing (e.g. fitting) separately, in case there is a problem. //DO YOUR STUFF HERE //CALL THIS LAST DSelector::Finalize(); //Saves results to the output file }