flyopunta/src/src_uti62/gigabeam.cpp

/***************************************************************************
 *   Copyright (C) 2005 by Giuseppe Pierazzini          Pisa 19.12.2010    *
 *   pierazzini@unipi.it                                                   *
 *                                     F l y o                             *
 *   MA62                                                                  *
 ***************************************************************************/
#include "parm.h"
#include "apparato.h"
#include "gvet.h"
#include "evento.h"
#include "reazione.h"
#include "clusgmp.h"
#include "gigabeam.h"
using namespace GmpConst;
using namespace std;

GigaBeam::GigaBeam()
{
    if (giga1==0||giga2==0||giga3==0)
    {
        Gout<<"\n ERROR....  Giga pointer not defined ---> exit(0)"<<endl;
        cout<<"\n ERROR....  Giga pointer not defined ---> exit(0)"<<endl;
        exit(0);

    }

    if (mgb1==0||mgb2==0||mgb3==0||mgb4==0||mgb5==0)

    {
        Gout<<"\n ERROR....  Giga pointer of magnet not defined ---> exit(0)"<<endl;
        cout<<"\n ERROR....  Giga pointer of magnet not defined ---> exit(0)"<<endl;
        exit(0);
    }




//   double ZG21,ZG32,ZG31,ZM21,ZG3M5;


    ZG21=giga2->Get_Cface().z-giga1->Get_Cface().z;
    ZG32=giga3->Get_Cface().z-giga2->Get_Cface().z;
    ZG31=giga3->Get_Cface().z-giga1->Get_Cface().z;
    ZM21 = mgb2->Get_Cface().z-mgb1->Get_Cface().z;  //==  ZM43
    ZG3M5=giga3->Get_Cface().z-mgb5->Get_Cface().z;
    Magkik1=mgb1->Get_Kik();  // lo stesso in valore assoluto per mgb2,mgb3,mgb4
    Magkik5=mgb5->Get_Kik();

};

int GigaBeam::GetCoordinate ()
{
    // g1,g2,g3  tutti >=0;
    // se pipeline è on si  verificano i tempi degli hits
    // se non corrispondono al tempo del trigger si  ritorna con errore -1
    if ( giga1->mhit<1|| giga2->mhit<1||giga3->mhit<1) return -1 ;
    
    

    int g1,g2,g3;
    g1=g2=g3=0;

// test on pipeline....if not skip to the end with g_index zero
    if (Pipeline==1)
    {
        int ok=0;
        double Spread=  giga1->Get_Res()*3.;  // stesso valore  per giga2 e giga3
        double t_ev = evento_.Gen.tempo;

        for ( g1=0;g1<giga1->mhit;g1++)
        {
            for (  g2=0;g2<giga2->mhit;g2++)
            {
                for ( g3=0;g3<giga3->mhit;g3++)
                {

                    if (fabs(giga1->M_Hits[g1].tempo-t_ev) > Spread  ||
                            fabs(giga2->M_Hits[g2].tempo-t_ev) > Spread ||
                            fabs(giga3->M_Hits[g3].tempo-t_ev) > Spread       )
                        continue;
                    //else
                    ok=1;
                    break;

                }
            }
        }
        if (ok==0) return -1; // no good hits....

    }

    xg1=giga1->M_Hits[g1].Xlab.x;
    xg2=giga2->M_Hits[g2].Xlab.x;
    xg3=giga3->M_Hits[g3].Xlab.x;

    yg1=giga1->M_Hits[g1].Xlab.y;
    yg2=giga2->M_Hits[g2].Xlab.y;
    yg3=giga3->M_Hits[g3].Xlab.y;

    return 0;


};
gvet  GigaBeam::Beam_Rec()
{

    gvet PBeam;
    PBeam.setvn(-0.09,0.,75.);
    
// Attenzione di solito si prende sempre il primo hits... nel caso di pipeline on
// si fa  una loop sugli hits disponibili: Vedi GetCoordinate()
// GetCoordinate fa un test sulla compatibilità temporale e accetta
// la tripletta che sta nella risoluzione temporale del Gigatrack
// altrimenti ritorna con -1

    if (GetCoordinate ()<0)   return PBeam;
    //OK
    // calcolo l'impuso ppz
    double dydz=(yg2-yg1)/ZG21-(yg3-yg2)/ZG32;
    double dzmdz= ZM21*(1./ZG21+1./ZG32);
    double ppz= Magkik1 *dzmdz/dydz;

    //calcolo di ppy
    double ppy = ppz*(yg3-yg1)/ZG31;
    //calcolo di ppx
    double ppx= ppz*(xg3-xg1)/ZG31 - Magkik5*ZG3M5/ZG31;
    ppx+=Magkik5;
    PBeam.setvn(ppx,ppy,ppz);
    return PBeam;
};