---

AlToolBox.cpp

Go to the documentation of this file.

//
// Implementation of the Tool Box.
//
// Author :  C. Delaere
//

#include "AlToolBox.h"
#include "AlephCollection.h"
#include "fortran_def.h"
#include "AlBjet.h"
#include "AlEflw.h"
#include "AlTrack.h"
#include "AlphaBanks.h"
#include "AlJet.h"

pair<AlephCollection<AlBjet>,vector<float> > AlToolBox::TagBjets(AlphaBanks* banks )
{
        int Error =0;
        int NumberOfJets =0;
        int NumberOfTracks =0;
        float Xjets[40] ={0};
        float Yjets[40] ={0};
        float Zjets[40] ={0};
        float Ejets[40] ={0};
        int TracksUsed[400] ={0};
        float TracksProba[400] ={0};
        int EflwsInJets[400] ={0};
        int JetsOfEflws[400] ={0};
        float JetProba[40] ={0};
        float HemiProba[2] ={0};
        float EventProba =0;
        F_QIPBTAG(Error,NumberOfJets,NumberOfTracks,Xjets,Yjets,Zjets,Ejets,TracksUsed,
                        TracksProba,EflwsInJets,JetsOfEflws,JetProba,HemiProba,EventProba);
        pair<AlephCollection<AlBjet>,vector<float> > output;
        AlephCollection<AlBjet> thejets;
        vector<float> values(banks->TrackPV().size()+1);
        if(Error==0)
        {
                values[values.size()-1] = EventProba;
                for(int i=0;i<NumberOfTracks;i++)
                {
                        if (TracksUsed[i]<=(int)values.size())
                                values[TracksUsed[i]-1] = TracksProba[i];
                        else    values.push_back(TracksProba[i]);
                }
                for(int i=0;i<NumberOfJets;i++)
                {
                        AlBjet myBjet;
                        myBjet.setA4V(Xjets[i],Yjets[i],Zjets[i],Ejets[i]);
                        myBjet.setBprobability(JetProba[i]);
                        myBjet.setScheme(1);
                        myBjet.setMetric(0);
                        thejets.push_back(myBjet);
                }
                for(unsigned int i=0;i<banks->EflwPV().size();i++)
                        if (JetsOfEflws[i]!=0) 
                                thejets[JetsOfEflws[i]-1].addObject(banks->EflwPV()[EflwsInJets[i]]);
        }
        output.first=thejets;
        output.second=values;
        return output;
}
AlephCollection<AlBjet> AlToolBox::ComputeBprobability(AlephCollection<AlJet*>* thejets, AlphaBanks* banks,bool negprob=false)
{
        // call QIPBTAG to get tracks probabilities
        pair<AlephCollection<AlBjet>,vector<float> > qipbtags = TagBjets( banks );
        // computes log factorials
        float logfac[401] = {0};
        for(int i=1;i<=400;i++) logfac[i] = logfac[i-1] + log(i);
        // some initialisations and declarations
        float sumjet[thejets->size()];
        int   ntrk_jet[thejets->size()];
        for (unsigned int i=0;i<thejets->size();i++) sumjet[i]=ntrk_jet[i]=0;
        // computes de probabilities products
        AlephCollection<AlTrack*> allthetracks = banks->TrackPV();
        AlephCollection<AlEflw*>  alltheeflws  = banks->EflwPV();
        int indexofthetrack = 0;
        for( AlephCollection<AlTrack*>::iterator itrack = allthetracks.begin(); itrack<allthetracks.end();itrack++,indexofthetrack++)
        {
                AlEflw* theasseflw = (*itrack)->getEflw();
                int ijet = 0;
                for(AlephCollection<AlJet*>::iterator jet = thejets->begin(); jet<thejets->end();jet++,ijet++)
                {
                        vector<AlObject*>  thejeteflws = (*jet)->getObjects();
                        for(vector<AlObject*>::iterator ieflw = thejeteflws.begin(); ieflw<thejeteflws.end();ieflw++)
                                if (((AlEflw*)(*ieflw))==theasseflw)
                                {
                                        float thetrack = qipbtags.second[indexofthetrack];      
                                        if ((negprob)&&(thetrack!=0))
                                        {
                                                sumjet[ijet] += ((thetrack>0) ? log(thetrack/2) : log(1+thetrack/2));
                                                ntrk_jet[ijet]++;
                                        }
                                        else if (thetrack>0)
                                        {
                                                sumjet[ijet] += log(thetrack);
                                                ntrk_jet[ijet]++;
                                        }
                                }
                }
        }
        // normalize the probability product for the number of measurements made
        int ijet=0;
        float loginvlog;
        float proba = 0;
        AlephCollection<AlBjet> result;
        for(AlephCollection<AlJet*>::iterator jet = thejets->begin(); jet<thejets->end();jet++,ijet++)
        {
                if(sumjet[ijet]<0)
                {
                        if((ntrk_jet[ijet]-1)>=1)
                                loginvlog = log(-sumjet[ijet]);
                        for(int i=1;i<ntrk_jet[ijet];i++)
                                proba += exp(i*loginvlog-logfac[i]);
                        float logproba = log(proba) + sumjet[ijet];
                        proba = exp((logproba>-50) ? logproba : -50) ;
                        proba = ((proba<1) ? proba : 1) ;
                }
                else
                        proba = 1;
                result.push_back(AlBjet(**jet,proba));
        }
        // done
        return result;
}

---