flyopunta/src/src_uti48/flyline.cpp

/***************************************************************************
 *   Copyright (C) 2005 by Giuseppe Pierazzini          Pisa 20.02.05      *
 *   pierazzini@unipi.it                                                   * 
 *                                     F l y o                             *
 *   Epsi/NA48                                                             *
 ***************************************************************************/
#include <math.h>
#include "flyoh.h"

//using namespace std;

// track finder for one  spectrometer NA48 apparatus

Flyline::Flyline()
{
  // per il momento si pensa anche al nuovo esperimento??? kpivv
  Title= " Ricostruzione delle tracce: nuova procedura base";
  Gout<<"\n"<<Title;
  First_trk=0;
  Tot_trk=0;   //totale di particelle ricostruite
  Totrac=0;
  for(int i=0;i<5;i++)
    {
      index[i]=0;
      xct[i]=yct[i]=zct[i]=vxt[i]=vyt[i]=vzt[i]=0.0;
      chiq[i]=ptrk[i]=charge[i]=0.0;
    }

  mas_pip=ptr_Pip->Get_Massa();
  // campo magnetico supposto lungo y di segno Bsegno..
  Bsegno1=Bsegno2=1.;
  kick1=kick2=0.0;

  if(Mag48!=0)   // esperimenti apparato vecchio
    {
      if(Mag48->Get_Field()[1]<0.)
        Bsegno1=-1.;  //fix the magnetic  field segno respect to Y axis
      kick1=Mag48->Get_Kik();
      Gout<<"\n B  Mag48 "<<Mag48->Get_Field()[1]<<"   segno  "<<Bsegno1<<"  kick "<< kick1;
    }
  else if(mnp1!=0&&mnp2!=0)  // apparato nuovo a 6 camere
    {
      if(mnp1->Get_Field()[1]<0.)
        Bsegno1=-1.;  //fix the magnetic  field segno respect to Y axis
      kick1=mnp1->Get_Kik();
      if(mnp2->Get_Field()[1]<0.)
        Bsegno2=-1.;  //fix the magnetic  field segno respect to Y axis
      kick2=mnp2->Get_Kik();
      Gout<<"\n B  mnp1 "<<mnp1->Get_Field()[1]<<"   segno  "<<Bsegno1<<"  kick "<< kick1;
      Gout<<"\n B  mnp2 "<<mnp2->Get_Field()[1]<<"   segno  "<<Bsegno2<<"  kick "<< kick2;
    }
  else
    {
      Gout<<" Field and kik error ... seen flyline.cpp  ... exit"<<std::endl;
      exit(0);
    }

}
Flyline::~Flyline()
{}

void  Flyline::Get_Hits_from_Dev()
{
  // procedure to get hits from the chamber.
  // It is general for four or three chambers...
  // the hits in dev coordinates are stored in xhit vectors
  // Gout<<"\n Ev: "<< evento_.Event<<"  "<<"Estraggo le coordinate dalle camere...";
  int i;
  // remeber : no hit nhit=-1;
  nhit1=nhit2=nhit3=nhit4=nhit5=nhit6=-1;

  if(Cam1!=0)
    {
      nhit1=Cam1->nhit;
      for(i=0;i<nhit1+1;i++)
        xhit1[i]=Cam1->M_Hits[i].Xlab;
    }
  if(Cam2!=0)
    {
      nhit2=Cam2->nhit;
      for(i=0;i<nhit2+1;i++)
        xhit2[i]=Cam2->M_Hits[i].Xlab;
    }
  if(Cam3!=0)
    {
      nhit3=Cam3->nhit;
      for(i=0;i<nhit3+1;i++)
        xhit3[i]=Cam3->M_Hits[i].Xlab;
    }
  if(Cam4!=0)
    {
      nhit4=Cam4->nhit;
      for(i=0;i<nhit4+1;i++)
        xhit4[i]=Cam4->M_Hits[i].Xlab;
    }
 
 //  Gout<<"\n Ev "<<evento_.Event<<"  Cordinate raccolte..."<<
 //        nhit1<<" "<<nhit2<<" "<<nhit3<<" "<<nhit4<<" "<<nhit5<<" "<<nhit6;
}
void Flyline::Find_Tracce()
{
  Tracc[0]->Reset();
  Tracc[1]->Reset();
  Tracc[2]->Reset();
  Tracc[3]->Reset();
  Tracc[4]->Reset();
  Tracc[5]->Reset();



  // get out if some important chamber has not hits....
  Tot_trk=-1;
  i5=i6=0;
  Tot_Charge=0.;
  Totrac=0;
  if( nhit1<0||nhit2<0)
    return;
  if( nhit3<0||nhit4<0)
    return;

  for(int i=0;i<6;i++)
    {
      for(int j=0;j<10;j++)
        {
          used_trk[i][j]=0;
        }
    }

  // find the tracks... max 3
  for(int trk=0;trk<6;trk++)
    {
      Iter_trk=0;
      TrkNata=0;
      chiqlast=300.;   // Chiq trk cut : to be tuned.
      //  Gout<<"\n Ev "<<evento_.Event<<"  ===================    track finding  "<<trk+1<<
      //       "  === <"<<nhit1<<nhit2<<nhit3<<nhit4<<"> ===========";
      // loop to find the tracks...
      //for(i1=nhit1;i1>=-1;i1--)   // attenzione i1==-1 accepted to fake a lack of hit..
      // {if(i1>-1){if(used_trk[0][i1]>0)continue;    // puo' mancare for I1==nhit1

      for(i1=nhit1;i1>-1;i1--)
        {
          if(used_trk[0][i1]>0)
            continue;
          y1=xhit1[i1][1];
          z1=xhit1[i1][2];
          for(i2=nhit2;i2>-1;i2--)
            {
              if(used_trk[1][i2]>0)
                continue;         // non puo' mancare???
              y2=xhit2[i2][1];
              z2=xhit2[i2][2];
              for(i3=nhit3;i3>-1;i3--)
                {
                  if(used_trk[2][i3]>0)
                    continue;       // non puo' mancare in kpivv experimento
                  y3=xhit3[i3][1];
                  z3=xhit3[i3][2];
                  for(i4=nhit4;i4>-1;i4--)
                    {
                      if(used_trk[3][i4]>0)
                        continue;      // non puo' mancare
                      y4=xhit4[i4][1];
                      z4=xhit4[i4][2];

                      double mtg=(y4-y1)/(z4-z1);
                      //        double dy1 = y1-y1-mtg*(z1-z1);
                      double dy2 = y2-y1-mtg*(z2-z1);
                      double dy3 = y3-y1-mtg*(z3-z1);
                      //        double dy4 = y4-y1-mtg*(z4-z1);

                      chiq_ora=(dy2*dy2+dy3*dy3)*1700.;

                      

                      Iter_trk++;
                      if(chiq_ora<chiqlast)    //chiq good.. save hits...
                        { chiqlast=chiq_ora;
                          Save_Traccia();
                          TrkNata=1;
                          if(chiq_ora<15.)
                            {
                              TrkNata=10;
                              break;
                            }
                        }

                    }
                  if(TrkNata==10)
                    break;
                }
              if(TrkNata==10)
                break;
            }
          if(TrkNata==10)
            break;
        }
      // }}
      if( TrkNata==0)
        break;
      Tot_trk++;

      used_trk[0][isv1]=1;
      used_trk[1][isv2]=1;
      used_trk[2][isv3]=1;
      used_trk[3][isv4]=1;

      Tracc[Tot_trk]->Get_Traccia();
      if(Tot_trk>4)
        break;
    }
  Totrac=Tot_trk+1;
}

void Flyline::Save_Traccia()
{
  isv1=i1;
  isv2=i2;
  isv3=i3;
  isv4=i4;
  isv5=i5;
  isv6=i6;
  IterSav=Iter_trk;
  chiqsav=chiq_ora;
  indexsav=i6*100000+i5*10000+i4*1000+i3*100+i2*10+i1;

  //  Gout<<"\n Selected: "<<" "<<isv1<<" "<<isv2<<" "<<isv3<<" "<<isv4<< setprecision(4)<<fixed<<
  //    "  chiqsav "<<chiqsav<<endl;
}