ShowRec_Alg_AG.cpp File Reference

This graph shows which files directly or indirectly include this file:

Functions
void	ReconstructShowers_AG ()

Function Documentation

ReconstructShowers_AG()

void ReconstructShowers_AG

(

)

{
    Log(2, "ShowRec.cpp", "--- void ReconstructShowers_AG() ---");
 
    //-----------------------------------------------------------------
    // Main function for reconstruction of "Advanced Gamma" Algorithm
    //-----------------------------------------------------------------
 
    //-----------------------------------
    // For each InitiatorBT this is
    // divided in several small parts:
    //
    // 1) Make local_gAli with cut parameters, Make GetPID of InBT and corresponding of plates
    // 3) Loop over (whole) local_gAli, check BT for Cuts
    // 4) Calculate pur/eff/NBT numbers
    // 5) Fill Trees
    //-----------------------------------
 
    // Define Helper Variables:
    EdbPVRec* local_gAli;
    EdbSegP* InBT;
    EdbSegP* seg;
    Float_t local_gAli_pat_interim_halfsize=0;
 
    GLOBAL_InBTArrayEntries=GLOBAL_InBTArray->GetEntries();
    GLOBAL_INBTSHOWERNR=0;
 
    //-----------------------------------------------------------------
    // Since GLOBAL_InBTArray is filled in ascending ordering by zpositon
    // We use the descending loop to begin with BT with lowest z first.
    for (Int_t i=GLOBAL_InBTArrayEntries-1; i>=0; --i) {
 
        //-----------------------------------
        // CounterOutPut
        if (gEDBDEBUGLEVEL==2) if ((i%1)==0) cout << GLOBAL_InBTArrayEntries <<" InBT in total, still to do:"<<Form("%4d",i)<< "\r\r\r\r"<<flush;
        //-----------------------------------
 
        //-----------------------------------
        // Get InitiatorBT from GLOBAL_InBTArray
        InBT=(EdbSegP*)GLOBAL_InBTArray->At(i);
        //--------
        GLOBAL_InBT_E=InBT->P();
        GLOBAL_InBT_TanTheta=TMath::Sqrt(InBT->TX()*InBT->TX()+InBT->TY()*InBT->TY());
        GLOBAL_InBT_Flag=InBT->Flag();
        GLOBAL_InBT_MC=InBT->MCEvt();
        //--------
        Int_t local_NBT=0;
        Int_t local_NBTMC=0;
        Int_t local_NBTallMC=0;
        Int_t local_NBTeMC=0;
        float_t local_pure=-1;
        float_t local_purall=-1;
        Int_t npat_int=0;
        Int_t npat_total=0;
        Int_t npatN=0;
        Int_t npat_Neff=0;
        Int_t NBT_Neff=0;
        Int_t NBTMC_Neff=0;
        Int_t NBTMCe_Neff=0;
        //--------
 
        if (gEDBDEBUGLEVEL>2) {
            cout << endl << endl << "--- Starting Shower for Number " << i << " now: "<<endl;
            InBT->PrintNice();
        }
        //-----------------------------------
 
        //-----------------------------------
        // 1) Make local_gAli with cut parameters:
        //-----------------------------------
        local_gAli = TransformEdbPVRec(GLOBAL_gAli, InBT);
        // Add InBT to GLOBAL_ShowerSegArray
        GLOBAL_ShowerSegArray -> Add(InBT);
        //-----------------------------------
 
 
        //-----------------------------------
        // 2) Loop over (whole) local_gAli, check InitiatorBT
        // 2) compatible with a segment forming a e+e- pair:
        //-----------------------------------
        Int_t NPairs=0,PIDDIFF=0,minPIDDIFF=0;
 
        Int_t local_gAli_npat=local_gAli->Npatterns();
        Int_t btloop_cnt_N=0;
        for (Int_t patterloop_cnt=local_gAli_npat-1; patterloop_cnt>=0; --patterloop_cnt) {
            if (gEDBDEBUGLEVEL>3) cout << "--- --- Doing patterloop_cnt= " << patterloop_cnt << endl;
            btloop_cnt_N=local_gAli->GetPattern(patterloop_cnt)->GetN();
            if (local_gAli->GetPattern(patterloop_cnt)->Z()<InBT->Z()) continue;
            if (TMath::Abs(local_gAli->GetPattern(patterloop_cnt)->Z()-InBT->Z())>3000) continue;
            if (local_gAli->GetPattern(patterloop_cnt)->Z()<InBT->Z()) continue;
            for (Int_t btloop_cnt=0; btloop_cnt<btloop_cnt_N; ++btloop_cnt) {
                seg = (EdbSegP*)local_gAli->GetPattern(patterloop_cnt)->GetSegment(btloop_cnt);
                if (gEDBDEBUGLEVEL>3) seg->PrintNice();
                if (IsSameSegment(seg,InBT)) continue;
                PIDDIFF=TMath::Abs(seg->PID()-InBT->PID());
                if (GetMinimumDist(InBT,seg)>50) continue;
                if (GetdR(InBT,seg)>200) continue;
                if (GetdeltaThetaSingleAngles(InBT,seg)>0.5) continue;
                if (PIDDIFF>minPIDDIFF) continue;
                //     cout << " GetMinimumDist(InBT,seg); " << GetMinimumDist(InBT,seg) << "PIDDIFF "<< PIDDIFF<< endl;
 
                AddBTToArrayWithCeck(seg, GLOBAL_ShowerSegArray);
                minPIDDIFF=PIDDIFF;
            }
        }
        PrintShowerObjectArray(GLOBAL_ShowerSegArray);
        //   return;
 
        //-----------------------------------
        // 2) Loop over (whole) local_gAli, check BT for Cuts
        //-----------------------------------
        if (gEDBDEBUGLEVEL>2) cout << "--- local_gAli_npat=  " << local_gAli_npat << endl;
 
        // Loop over all plates of local_gAli, since this is already
        // extracted with the right numbers of plates...
        for (Int_t patterloop_cnt=local_gAli_npat-1; patterloop_cnt>=0; --patterloop_cnt) {
            if (gEDBDEBUGLEVEL>3) cout << "--- --- Doing patterloop_cnt= " << patterloop_cnt << endl;
            btloop_cnt_N=local_gAli->GetPattern(patterloop_cnt)->GetN();
 
            // if pattern is before the Initiator BT->Z() position we just go on:
            if (local_gAli->GetPattern(patterloop_cnt)->Z()<InBT->Z()) continue;
 
            for (Int_t btloop_cnt=0; btloop_cnt<btloop_cnt_N; ++btloop_cnt) {
                seg = (EdbSegP*)local_gAli->GetPattern(patterloop_cnt)->GetSegment(btloop_cnt);
                if (gEDBDEBUGLEVEL>3) seg->PrintNice();
 
                // Now apply cut conditions: Advanced Gamma Alg  --------------------
                if (!GetConeOrTubeDistanceToBTOfShowerArray(seg, InBT, GLOBAL_ShowerSegArray, CUT_PARAMETER[0], CUT_PARAMETER[1])) continue;
                if (!FindPrecedingBTsSingleThetaAngle(seg, InBT, local_gAli, GLOBAL_ShowerSegArray)) continue;
                // end of    cut conditions: Advanced Gamma Alg  --------------------
 
                // If we arrive here, Basetrack seg has passed criteria
                // and is then added to the shower array:
                // Check if its not the InBT which is already added:
                //     if (IsSameSegment(seg,InBT)) { ; } // do nothing;
                //     else {GLOBAL_ShowerSegArray -> Add(seg);}
                AddBTToArrayWithCeck(seg, GLOBAL_ShowerSegArray);
            }
 
 
            // Calc BT density around shower:
            EdbPattern* pat_interim=local_gAli->GetPattern(patterloop_cnt);
            CalcTrackDensity(pat_interim,local_gAli_pat_interim_halfsize,npat_int,npat_total,npatN);
 
            // Calc TrackNumbers for plate for efficency numbers:
            CalcEfficencyNumbers(pat_interim, InBT->MCEvt(), NBT_Neff, NBTMC_Neff,NBTMCe_Neff);
        }
        // end of loop over all plates of local_gAli
 
        if (gEDBDEBUGLEVEL>2) PrintShowerObjectArray(GLOBAL_ShowerSegArray);
 
 
 
 
        Int_t s_NBT=0;
        Int_t s_NBTMC=0;
        Int_t s_NBTallMC=0;
        Int_t s_NBTeMC=0;
        Double_t  s_eff=0;
        Double_t s_purall=0;
        Double_t s_pure=0;
        CalcEffPurOfShower(GLOBAL_ShowerSegArray, s_NBT, s_NBTMC, s_NBTallMC, s_NBTeMC, s_purall, s_pure);
 
        //-----------------------------------
        // 4) Calculate pur/eff/NBT numbers,
        // not needed when only reconstruction
        // done:
        //-----------------------------------
        if (cmd_OUTPUTLEVEL>=2 || cmd_OUTPUTLEVEL==0 ) {
            Int_t NBT=0;
            Int_t NBTMC=0;
            Int_t NBTallMC=0;
            Int_t NBTeMC=0;
            Double_t  eff, purall, pure;
            CalcEffPurOfShower2(GLOBAL_ShowerSegArray, NBT, NBTMC, NBTallMC, NBTeMC, purall, pure, NBT_Neff, NBTMC_Neff,NBTMCe_Neff);
 
            // Fill only for MC Event:
            if (GLOBAL_InBT_MC>0) {
                GLOBAL_EvtBT_Flag=GLOBAL_EvtBT_FlagArray[GLOBAL_InBT_MC];
                GLOBAL_EvtBT_MC=GLOBAL_EvtBT_MCArray[GLOBAL_InBT_MC];
                GLOBAL_EvtBT_E=GLOBAL_EvtBT_EArray[GLOBAL_InBT_MC];
                GLOBAL_EvtBT_TanTheta=GLOBAL_EvtBT_TanThetaArray[GLOBAL_InBT_MC];
            }
            GLOBAL_trckdens=shower_trackdensb;
        }
 
 
        //-----------------------------------
        // 5) Fill Tree:
        //-----------------------------------
        TREE_ShowRecEff->Fill();
        if (gEDBDEBUGLEVEL>3) TREE_ShowRecEff->Show(TREE_ShowRecEff->GetEntries()-1);
 
 
        //-----------------------------------
        // 6a) Transfer ShowerArray to treebranchTreeEntry:
        //-----------------------------------
        if (cmd_OUTPUTLEVEL>0) {
            TransferShowerObjectArrayIntoEntryOfTreebranchShowerTree(TREE_ShowShower,GLOBAL_ShowerSegArray);
        }
 
 
        //------------------------------------
        // Reset and delete important things:
        // also  to avoid memory problems ...
        //-----------------------------------
        GLOBAL_ShowerSegArray->Clear();
        if (gEDBDEBUGLEVEL>3) cout << "--- ---GLOBAL_ShowerSegArray->GetEntries(): "<< GLOBAL_ShowerSegArray->GetEntries() << endl;
        delete local_gAli;
        local_gAli=0;
        ++GLOBAL_INBTSHOWERNR;
        //------------------------------------
    }
    // end of loop over GLOBAL_InBTArrayEntries
    //-----------------------------------------------------------------
 
    if (gEDBDEBUGLEVEL==2) cout << endl<<flush;
    if (gEDBDEBUGLEVEL>3) cout << "---TREE_ShowRecEff->GetEntries() ... " << TREE_ShowRecEff->GetEntries() << endl;
    if (gEDBDEBUGLEVEL>3) cout << "---GLOBAL_INBTSHOWERNR ... " << GLOBAL_INBTSHOWERNR<< endl;
 
 
    return;
}