ShowRec_Alg_RC.cpp File Reference

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

Functions
void	ReconstructShowers_RC ()

Function Documentation

ReconstructShowers_RC()

void ReconstructShowers_RC

(

)

last to variables are dummyvariabes.

if this is in, it doesnt work for iteratet ShowerAxisCenterGravityBT

{
    Log(2, "ShowRec.cpp", "--- void ReconstructShowers_RC() ---");
 
    //-----------------------------------------------------------------
    // Main function for reconstruction of "Recursive Cone" 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 BT for Cuts
        //-----------------------------------
        Int_t local_gAli_npat=local_gAli->Npatterns();
        Int_t btloop_cnt_N=0;
        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: RECURSIVE CONE 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: RECURSIVE CONE 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:
                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 f_NBT=0;
        Int_t f_NBTMC=0;
        Int_t f_NBTallMC=0;
        Int_t f_NBTeMC=0;
        Double_t  f_eff=0;
        Double_t f_purall=0;
        Double_t f_pure=0;
        CalcEffPurOfShower(GLOBAL_ShowerSegArray, f_NBT, f_NBTMC, f_NBTallMC, f_NBTeMC, f_purall, f_pure);
 
        //-----------------------------------
        // 3) Try RE-ITERATION over (whole) local_gAli!!
        // 3) Using the iterated shower Gravity Center!
        // 3)
        // 3) WILL NOT WORK FOR ELECTRONS ... CAUSE FISRT POSITION IS SHIFTED AFTER ITERATION !!!
        // 3) DEPRECIAPTED
        // 3) DEPRECIAPTED  .... DO NOT USE IT NOW!!!!!!!!!!!!!!
        // 3) DEPRECIAPTED  (even if it works we dont need to reevaluate the
        // 3)                efficency numbers anymore, since calculated above)
        //-----------------------------------
 
        Bool_t DoIteration=kFALSE;
        //    Bool_t DoIteration=kTRUE// 3) DEPRECIAPTED  .... DO NOT USE IT NOW!!!!!!!!!!!!!!
 
        if (DoIteration) {
 
            cout << "-------------------------------------------"<<endl;
            EdbSegP* ShowerAxisCenterGravityBT = BuildShowerAxis(GLOBAL_ShowerSegArray);
            GLOBAL_ShowerSegArray->Clear();
            //-----------------------------------
            // 3) Loop again over (whole) local_gAli, check BT for Cuts
            //-----------------------------------
            local_gAli_npat=local_gAli->Npatterns();
            btloop_cnt_N=0;
            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();
 
                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: RECURSIVE CONE Alg  --------------------
                    //      if (!GetConeOrTubeDistanceToBTOfShowerArray(seg, ShowerAxisCenterGravityBT, GLOBAL_ShowerSegArray, CUT_PARAMETER[0], CUT_PARAMETER[1])) continue;
                    if (!FindPrecedingBTsSingleThetaAngle(seg, ShowerAxisCenterGravityBT, local_gAli, GLOBAL_ShowerSegArray)) continue;
                    // end of    cut conditions: RECURSIVE CONE 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 (seg->X()==ShowerAxisCenterGravityBT->X()&&seg->Y()==ShowerAxisCenterGravityBT->Y()) { ; } // do nothing;
                    //      else {GLOBAL_ShowerSegArray -> Add(seg);}
                    AddBTToArrayWithCeck(seg, GLOBAL_ShowerSegArray);
                }
            }
            // end of loop over all plates of local_gAli
            if (gEDBDEBUGLEVEL>2) PrintShowerObjectArray(GLOBAL_ShowerSegArray);
            cout << "-------------------------------------------"<<endl;
 
 
        } // of if (DoIteration)
 
        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);
        }
 
 
        //-----------------------------------
        // 7) cout some numbers from shower reco:
        //-----------------------------------
        if (gEDBDEBUGLEVEL>2) {
            cout << "-------------------------------  _NUMBERS ----------------------"<<endl;
            cout << "ITERATION 0 : --- NBT    = " << f_NBT << endl;
            cout << "ITERATION 0 : --- NBTMC  = " << f_NBTMC << endl;
            cout << "ITERATION 0 : --- purall = " << f_purall << endl;
            cout << "------------"<<endl;
            cout << "ITERATION 1 : --- NBT    = " << s_NBT << endl;
            cout << "ITERATION 1 : --- NBTMC  = " << s_NBTMC << endl;
            cout << "ITERATION 1 : --- purall = " << s_purall << endl;
            cout << "NUMBERS  " << f_NBT << " " <<  f_NBTMC << " " <<  f_purall << " " <<  s_NBT << " " <<  s_NBTMC << " " <<  s_purall << " " << endl;
        }
 
        //------------------------------------
        // 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;
}