|
FlyoDoc_2011 Pisa 2011 by GmP --- 011
|
|
FlyoDoc_2011 Pisa 2011 by GmP --- 011
|
#include <devflt.h>
Public Member Functions | |
| DevFlt () | |
| double | Get_CamMag () |
| virtual int | Traccia (Particella *) |
| ================ T R A C C I A ==================================== | |
| virtual void | Prgeom () |
Protected Attributes | |
| gvet | B |
| gvet | X |
| gvet | XIN |
| gvet | PIN |
| gvet | VSIN |
| gvet | Puls |
| gvet | KIK |
| gvet | Pvers |
| gvet | Dvers |
| int | iter_mag |
| double | step |
| double | dse |
| double | path |
| double | rho |
| double | freestep |
| double | vol_max |
| double | cam_max |
| DevFlt::DevFlt | ( | ) |
magnete Mag.con buco per filtro DevFlt 22 attenzione rivedre bene la definizione dei vettori dev....per rotazioni e traslazioni dal lab.....
Definition at line 27 of file devflt.cpp.
: DevRtCl() { devtype=TypDevFlt; devclass="Magn_Filter"; // campo costante campo=0 default; // campo variabile campo>0; Definire.... ricorda! campo=0; }
| int DevFlt::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 59 of file devflt.cpp.
{
// see the DevMagn:: Traccia()
// 18/03/03 update: general tracing.. the magnetic field now could be variable
// and without any restriction..
// attention! inside the magnet no brems, nor pairprod for the moment..
// the track tracing is iterative with step size small enough to keep the
// tracing errors to very low values
//
//*****
// La carica elementare nel sistema Ibrido e'
// ec=.299792458[ Gev/c * T**(-1) m **(-1) ]
// va qui moltiplicata per 10e-6 poiche'
// si usa il Gauss ed il Centimetro; segue...
// r(cm)=[p cos(lambda)(GeV)]/(2.9979258e-7 * H(Gauss) )
const double ec=2.99792458e-7;
//****
// qui il campo è differente da zero solo nel ferro... il campo nel buco è nullo!!!
// si vuole filtrare via i muoni....
if(nowpos==ESTERNO || nowpos==INBUCO) // se esterno o nel buco traccia regolarmente...
{ Device::Traccia(pr); return 0;}
// allora e' INTERNO: si traccia in campo magnetico.
if (pr->P.norma == 0.0) return 0; // no motion ....
if (pr->Get_Charg() == 0.0) // neutra traccio regolamente
{ Device::Traccia(pr); return 0;}
//------------
X=Lab2cDev(pr->X); // Coordinata di impatto sul dev
PIN=Puls=Lab2Dev(pr->P); // impulso vettore in device....
VSIN=Pvers=Puls.Verso();
XIN=X; // device frame coordinates
vol_max=pr->pathok; // cammino potenziale prima di decadere!
// P_Cin=X-Centrin; //Occorre in Get_CamMag(), specie per buchi non simmetrici
P_Cin=X;
freestep=cam_max=Get_CamMag(); // cammino potenziale in linea retta per uscire dal magnete
if(freestep>vol_max) freestep=vol_max;
/* Campo in Field
campo=0 B field costante (va definito nella classe)
campo>0 B field variabile
*/
if(campo==0)
{
if (Field.norma == 0.0)
{ Gout<<"\n *** error*** Standard magnet "<<nome<<" null field???";
exit(0);}
}
// se variabile ridefinisco il campo iniziale
else Get_Field(P_Cin);
// nota il sistema di riferimento dei vettori e' quello del dev.
// predico il raggio di curvatura.... massimo...
rho= 9999999. ;
if(Field.norma>0.0001)rho = Puls.Norma()/(ec*Field.norma);
// il calcolo corretto e' rho= (P&H).Norma/(H_norma**2*ec)
/*************/
if(rho<2.*freestep) return 1; // la particella resta intrappolata....
// definizione dello step.
step =rho/3000.;
if(step<2.)step=2.;
else if(step>25.)step=25.;
if(step>freestep) step= freestep+0.04;
dse = step*ec*pr->Get_Charg();
path=0.0;
iter_mag=0;
if(Debugon==1)
{// ok nel caso di campo costante ==> valore significativo
double forkik=ec*Field.norma*cam_max*pr->Get_Charg();
Gout<<"\n\n Ev="<<evento_.Gen.Event<<" "<<pr->Get_Nome()<<" in "<<nome
<<" pathok="<<pr->pathok<<" Cammx="<< cam_max
<<" freestep "<<freestep
<<" Kik="<<forkik<<" rho="<<rho;
Field.print("Field");
}
// tracing main loop =======================
while (iter_mag++<MAXITER)
{
if(campo>0) Get_Field(P_Cin); //variabile
/*
evento_.Event,nome,pr->Get_Nome(),pr->pathok,cam_max,rho);
Pvers.print("pvers");
P_Cin.print(pr->Get_Nome());
Field.print("Field");
*/
Puls+= Pvers&Field * dse; // Aggiorno P dopo un cammino ds
Pvers=Puls.NVerso(); // calcola la nuova direzione di volo
P_Cin+= step*Pvers; // aggiorno la posizione
path+=step; // aggiorno il cammino
// check for end iteration...
if(step>freestep) break; // decaduta...o uscita dal mag...
freestep-=step; //aggiorno il cammino da fare
if (fabs(X[0])> Lout[0] || fabs(X[1])> Lout[1] ||fabs(X[2])> Lout[2]) break;
// if (fabs(P_Cin[0])> Loin[0] || fabs(P_Cin[1])> Lin[1] ||fabs(P_Cin[2])> Lin[2]) break;
if(step>freestep) // se supero con uno step il cammino restante allora...
{
// allora calcolo il cammino potenziale restante...
freestep=Get_CamMag();
if(freestep>vol_max-path)freestep=vol_max-path;
if(step>freestep)
{step= freestep+0.04;
if(step<0.06) break;
dse = step*ec*pr->Get_Charg();
}
}
}
// tracing end ====================================
if (iter_mag> 499)
{
Gout<<"\n\n Ev="<<evento_.Gen.Event<<" Trac. in "<<nome
<<" iter "<<iter_mag<<" for "<< pr->Get_Nome()<< " P norma "<< pr->P.norma
<<" camx "<< cam_max<<" rho="<<rho<<" step "<<step;
return (1);
}
// ===> go back to the apparatus system
pr->P.setvn(Dev2Lab(Puls),pr->Get_Massa()); // nota calcola anche la norma
pr->Vers=pr->P.Verso();
X=P_Cin+Centrin;
pr->X=Devc2Lab(X)+ 0.01*pr->Vers; // non norma... + piccolo spostamento
pr->pathok=path;
pr->path_done+=path;
/*****************************/
if(Debugon==1)
{
KIK=Puls-PIN;
Devc2Lab(XIN).print("Xin");
pr->X.print("Xout");
PIN.print("Pin");
Puls.print("Pout");
KIK.Norma();
KIK.print("KIK");
// Field.print("Last Field");
Gout<<"\n.. "<< nome <<" for "<<pr->Get_Nome()<<" Patok="<<freestep<<" Path "
<<path<<" Iter="<< iter_mag<<std::endl;
}
/*****************************/
return (0);
}
