FlyoDoc_2011 Pisa 2011 by GmP --- 011
|
Public Member Functions | |
virtual int | Traccia (Particella *) |
================ T R A C C I A ==================================== | |
virtual void | DataSmear () |
virtual void | Prgeom () |
virtual int | SimulaDev () |
virtual void | Reset_RootData () |
int | Get_DeadCel (double x, double y) |
void | Pad_Filling () |
double | E_pio_Relesed (Bufdev *) |
Public Attributes | |
double | sepm |
int | padhit |
float | dist_max |
float | epad28 |
float | m28x |
float | m28qx |
float | m28y |
float | m28qy |
float | epad44 |
float | m44x |
float | m44qx |
float | m44y |
float | m44qy |
float | quad1 |
float | quad2 |
float | quad3 |
float | quad4 |
Protected Attributes | |
double | x |
double | y |
double | ep |
double | se |
double | sx |
int | i |
int | ht |
int | hi |
int | dead |
double | xyd |
double | ee |
double | rr |
double | xfr |
double | slope |
int | Fusione |
Static Protected Attributes | |
static const float | Deadcl [90][2] |
void DevLkr::DataSmear | ( | ) | [virtual] |
Ridefinisce le corrdinate e le energie rivelate tenendo conto degli errori. Per il calcolo dello smearing usa le formula aggiornate standard usate in NA31 per il calorimetro.
Reimplemented from Device.
Definition at line 492 of file devlkr.cpp.
{ // operate on the dev_variables of type (2) M_Hits // double x,y,ep,se,sx; const char *partnom; i=0; while ( i<=nhit ) { x=M_Hits[i].Xdev.x; y=M_Hits[i].Xdev.y; partnom=M_Hits[i].pnome; ep=M_Hits[i].e_rivela; //energia rilasciate nel dev // if ( partnom==Gam|| partnom==Elec|| partnom==Elep ) { //================================================================ // new model from eta runs fit // double rr, xfr=0.,slope=0.0; int Tail_todo=1; if ( Tail_todo==1 ) { // if ( Pran() <=0.085 ) // first small correction..to be checked if ( Pran() <=0.090 ) // 19.07.09 { if ( ep<=30.0 ) slope=32.5+1.25*ep; else slope=75.0; rr =Pran(); if ( rr>1.e-9 ) // if ( rr>1.e-9 ) // 19.07.09 { ep*= ( 1.0+log ( rr ) /slope ); SelRea->Tailsrt_done++; } } // add long tails from E/P nel 3 per mille if ( Pran() <=0.0021 ) // hadronproduction correction.. // if ( Pran() <=0.003 ) // 19.07.09 { rr =Pran(); if ( rr>=0.000245 ) // if ( rr>=0.0001 ) // 19.07.09 xfr=0.92+log ( rr ) /10.14; else xfr=.1; ep*=xfr; SelRea->Tailong_done++; } // tails done } //================================================================ // calcolo degli errori per il calorimetro // x=Xdev[i][0]; y=Xdev[i][1]; // se= sqrt(0.01 + (0.001024+2.5e-5*ep)*ep); // 1998 se= sqrt ( 0.0081 + ( 0.001024+1.764e-5*ep ) *ep ); //1999 ep=Pgauss ( se,ep ); sx=sqrt ( 3.6e-3 + 0.1764/ep ); // correzione 2001 // sx=sqrt ( 2.5e-3 + 0.1806/ep ); // 19.07.09 ... vedere nella ricostruzione.... // x=Pgauss ( sx,x ); y=Pgauss ( sx,y ); } else { // le altre particelle...pioni o "fuse" ..si usa un errore di 1 cm // per i mu si usa 2 mm. double sigma=1.; // for handron or "Fuse" if ( partnom==Mum||partnom==Mup ) sigma=0.2; x=Pgauss ( sigma,x ); y=Pgauss ( sigma,y ); se= sqrt ( 0.0081 + ( 0.001024+1.764e-5*ep ) *ep ); //1999 ep=Pgauss ( se,ep ); } M_Hits[i].e_rivela=ep; M_Hits[i].Xdev.x=x; M_Hits[i].Xdev.y=y; M_Hits[i].Xlab=Devc2Lab ( M_Hits[i].Xdev ); i++; } }
double DevLkr::E_pio_Relesed | ( | Bufdev * | MH | ) |
Procedura per il calcolo dell'energia rilasciata dal pi carico in funzione della penetrazione.
Definition at line 838 of file devlkr.cpp.
{ // nota MH sono i parametri misurati di un hit (vedi BufDev *) static double anor=0.0; // define the energy left in the hadronic interction if it has flagged idead=-3 if ( MH->idead!=-3 ) return 0.0; double brnd=Pran(); if ( brnd<.3 ) return 0.; double rp=0.0, fp=0.0; //double alfa=10., beta=10.; originali double alfa=12., beta=10.; if ( anor==0.0 ) // prima volta calcolo l'integrale della fun. di distribuzione. { while ( rp<1. ) { anor+= ( 1.-exp ( -alfa*rp ) ) *exp ( -beta*rp ) *rp*rp; rp+=0.01; } } // normal path(gvet_0) double rnd=Pran(); rp=0.0; while ( rp<1. ) { fp+= ( 1.-exp ( -alfa*rp ) ) *exp ( -beta*rp ) *rp*rp/anor; if ( fp>rnd ) break; rp+=0.01; } // correzione dell'energia rilasciata in funzione della profondita' di interazione double depth= ( 2.*Lout.z-MH->last_path ) /60.; double reduci =1.- exp ( -depth*depth ); return rp* MH->Pdev.e*reduci; }
void DevLkr::Prgeom | ( | ) | [virtual] |
Stampa dei parametri del calorimetro.
Reimplemented from DevCrn.
Definition at line 483 of file devlkr.cpp.
{ Gout<<"================== R i v e l a t o r e ================="; Gout<<"\n --> < "<<nome <<" > [ "<<fun <<" ] Calorim. "; Device::Prgeom(); }
int DevLkr::SimulaDev | ( | ) | [virtual] |
Qui si simula il calorimetro per rendere le misure più realistiche. Per esempio due gamma troppi vicini (2 cm) si fondono in un solo gamma, o un gamma viene eliminato se cade su di una cella morta.
=================
Reimplemented from Device.
Definition at line 611 of file devlkr.cpp.
{ // this procedure performs a filter of the measures // to simulate the device response // attention for no gamma particles... // opera sulla variabili di tipo (2) predifinite M_Hits if ( mhit<1 ) return 0; emDev=0.; dead=0; // attenzione le dimensioni sono legate al massimo numere dellle particelle viste dal dev.... int flag[30]={0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}; // hit number of gammas in the associated det. (remember:-1, 0,1,2,3,4,5) const char *pnome; // ---------- // set flags and for pions try to make life a little bit worse into lkr i=0; while ( i<=nhit ) { flag[i]=1; // set all flags to 1 pnome=M_Hits[i].pnome; if ( pnome==Pip||pnome==Pim ) { M_Hits[i].e_rivela+=E_pio_Relesed( &M_Hits[i] ); // released pion energy } i++; } // cuts double Ecutg=0.35; // double Ecutg=.35; ht=nhit; while ( ht>-1 ) { // qui tagli in energia sui gamma o elettroni pnome=M_Hits[ht].pnome; if ( ! ( pnome==Gam||pnome==Elec||pnome==Elep ) ) { ht--; continue; } if ( flag[ht]!=0&& ( M_Hits[ht].e_rivela<Ecutg ) ) { flag[ht]=-2; dead++; } // geometry else if ( ( xyd= M_Hits[ht].Xdev.XYNorma() ) < 10. ) //taglio sul tubo del calorimetro { flag[ht]=-1; dead++; } // geometrical cut: else //if still alive, then dead cels test.... { x= M_Hits[ht].Xdev.x; y= M_Hits[ht].Xdev.y; int killcel=Get_DeadCel ( x,y ); if ( killcel==1 ) { flag[ht]=-3; dead++; } } ht--; } // cout<<"\n Lkr nhit dead "<< nhit<<" "<<dead<<" "<<killcel<<" "<<flag[0]; // for (int i=0;i<=nhit;i++) // Gout<<"\n Lk: "<< evento_.Gen.Event<<setprecision(3) // <<" xy "<< setw(8)<<M_Hits[i].Xdev.x<<" "<< setw(8)<<M_Hits[i].Xdev.y // <<" e "<< setw(8)<< M_Hits[i].e_rivela<<" pr "<<M_Hits[i].idm<< " flg "<<flag[i]; // Fusion procedure. Default to 1 ... zero for background test.... if ( Fusione==1 ) { evento_.Gen.Ifus=0; double ee1=0.0,ee2=0.0,eet=0.0; ht=nhit; while ( ht>-1 ) { if ( flag[ht]<1 ) { ht--; continue; } ee1=M_Hits[ht].e_rivela; gvet Daux=M_Hits[ht].Xdev; // Daux.print("Daux"); // hi=ht-1; while ( hi>-1 ) { if ( flag[hi]<1 ) { hi--; continue; } //skip... sepm = ( M_Hits[hi].Xdev-Daux ).XYNorma(); // distanza proiettata su xy // if ( sepm<10. ) return 0; // attenzione per Kl->pio0 ???? if ( sepm<2. ) // allora fusion sepm<2. { ee2=M_Hits[hi].e_rivela; eet=ee1+ee2; if ( eet<0.00001 ) { M_Hits[hi].Xdev= ( M_Hits[ht].Xdev+ M_Hits[hi].Xdev ) *0.5; M_Hits[hi].Xlab= ( M_Hits[ht].Xlab+ M_Hits[hi].Xlab ) *0.5; M_Hits[hi].tempo= ( M_Hits[ht].tempo + M_Hits[hi].tempo ) *0.5; } else { M_Hits[hi].Xdev= ( M_Hits[ht].Xdev*ee1+ M_Hits[hi].Xdev*ee2 ) /eet; M_Hits[hi].Xlab= ( M_Hits[ht].Xlab*ee1+ M_Hits[hi].Xlab*ee2 ) /eet; M_Hits[hi].tempo= ( M_Hits[ht].tempo*ee1 + M_Hits[hi].tempo*ee2 ) /eet; } M_Hits[hi].Pdev+=M_Hits[ht].Pdev; M_Hits[hi].Vers=M_Hits[hi].Pdev.Verso(); M_Hits[hi].e_rivela=eet; // to remeber that is an artefact define idm =-1 M_Hits[hi].idm=-1; if ( evento_.Gen.Event<1000 ) Gout<<" \n Fuse: "<<setw ( 8 ) <<evento_.Gen.Event<<setprecision ( 3 ) <<" "<<setw ( 8 ) << M_Hits[hi].pnome<<" di "<<setw ( 8 ) <<ee2 <<" con " <<M_Hits[ht].pnome <<" di "<<setw ( 8 ) <<ee1<<" sp "<<setw ( 8 ) <<sepm<<" " << M_Hits[hi].idm; M_Hits[hi].pnome=Fuse; flag[ht]=-4; dead++; evento_.Gen.Ifus++; break; } hi--; } ht--; } } // end fusione // soglia rivelazione nel kripton 400. MeV for ( int i=0;i<=nhit;i++ ) { if ( flag[i]>=0 ) { if ( M_Hits[i].e_rivela< 0.35 ) { flag[i]=-5; dead++; } else { double xcl=M_Hits[i].Xdev.x; double ycl=M_Hits[i].Xdev.y; if ( ( xcl*xcl+ycl*ycl ) <144. ) { flag[i]=-6; dead++; } } } } // compact the mesures.... int htt=-1; for ( int i=0;i<=nhit;i++ ) { if ( flag[i]>=0 ) { htt++; emDev+=M_Hits[htt].e_rivela= M_Hits[i].e_rivela; if ( htt==i ) continue; //avoid to rewrite the same data into itsself... M_Hits[htt]=M_Hits[i]; } } mhit=htt+1; //attenzione ridefinizione da 1 a n anziche' da 0 a n-1;/* nhit=mhit-1; //cout<<"\n Lkr mhit "<< nhit<<" "<<mhit; // pad filling per trigger zero level Pad_Filling( ); //================================= if ( Debugon==1 ) Gout<<"\n CalSimul "<<evento_.Gen.Event <<" mhit "<<setw ( 3 ) <<mhit<<" of "<<nhit+1 <<" dead "<<setw ( 3 ) << dead <<" fusi "<<setw ( 3 ) << evento_.Gen.Ifus <<" flag "<< flag[0]<<" "<<flag[1]<<" "<<flag[2] <<" " << flag[3]<<" "<<flag[4]<<" "<<flag[5]<<endl; return htt; }
int DevLkr::Traccia | ( | Particella * | pr | ) | [virtual] |
================ T R A C C I A ====================================
Traccia le particelle nell'apparato per un tratto come definito in Get_camm(). Nel tracciare la particella si tiene conto della possibilita' della pair production e del bremsstrhalung.
Reimplemented from Device.
Definition at line 449 of file devlkr.cpp.
{ // controllo il tracing della particella // se è un pione if ( pr->Get_Nome() ==Pip||pr->Get_Nome() ==Pim ) { // ed è interno al lkr ed ha un cammino potenziale maggiore di un cm.. if ( nowpos==INTERNO&&pr->pathok>1. ) { // calcolo in cammino prima di interagire double pathfai = -liint* log ( EPSR + Pran() ); if ( pathfai<pr->pathok ) { // interagisce prima di fare il cammino potenziale // restante dichiaro la particella come distrutta: Ifato=idead=-3; pr->Set_Fato ( Distr,-3 ); pr->pathok=pathfai; // ridefinisco il cammino potenziale nel dev prima dellinterazione distruttiva... } } } // se non è un pione traccio seguendo la procedura standard... // se è un pione traccio la parte di cammino restante prima della distruzione.... int ret=Device::Traccia ( pr ); return ret; }