12 |
extern "C" { |
extern "C" { |
13 |
void dotrack_(int*, double*, double*, double*, double*, int*); |
void dotrack_(int*, double*, double*, double*, double*, int*); |
14 |
void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*); |
void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*); |
15 |
int readb_(const char*); |
// int readb_(const char*); |
16 |
|
int readb_(); |
17 |
void mini2_(int*,int*,int*); |
void mini2_(int*,int*,int*); |
18 |
void guess_(); |
void guess_(); |
19 |
} |
} |
21 |
// |
// |
22 |
// |
// |
23 |
//-------------------------------------- |
//-------------------------------------- |
24 |
|
/** |
25 |
|
* Evaluates the trajectory in the apparatus associated to the track. |
26 |
|
* It integrates the equations of motion in the magnetic field. The magnetic field should be previously loaded ( by calling TrkLevel2::LoadField() ), otherwise an error message is returned. |
27 |
|
* @param t pointer to an object of the class Trajectory, |
28 |
|
* which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ). |
29 |
|
* @return error flag. |
30 |
|
*/ |
31 |
|
int Trajectory::DoTrack2(float* al){ |
32 |
|
|
33 |
|
double *dxout = new double[t->npoint]; |
34 |
|
double *dyout = new double[t->npoint]; |
35 |
|
double *dthxout = new double[t->npoint]; |
36 |
|
double *dthyout = new double[t->npoint]; |
37 |
|
double *dtlout = new double[t->npoint]; |
38 |
|
double *dzin = new double[t->npoint]; |
39 |
|
double dal[5]; |
40 |
|
|
41 |
|
int ifail = 0; |
42 |
|
|
43 |
|
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
44 |
|
for (int i=0; i<t->npoint; i++) dzin[i] = (double)z[i]; |
45 |
|
|
46 |
|
dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
47 |
|
|
48 |
|
for (int i=0; i<t->npoint; i++){ |
49 |
|
x[i] = (float)*dxout++; |
50 |
|
y[i] = (float)*dyout++; |
51 |
|
thx[i] = (float)*dthxout++; |
52 |
|
thy[i] = (float)*dthyout++; |
53 |
|
tl[i] = (float)*dtlout++; |
54 |
|
} |
55 |
|
|
56 |
|
return ifail; |
57 |
|
}; |
58 |
|
//--------------------------------------------- |
59 |
|
//--------------------------------------------- |
60 |
TrkTrack::TrkTrack(){ |
TrkTrack::TrkTrack(){ |
61 |
seqno = -1; |
seqno = -1; |
62 |
image = -1; |
image = -1; |
286 |
/** |
/** |
287 |
* Load the magnetic field |
* Load the magnetic field |
288 |
*/ |
*/ |
289 |
void TrkTrack::LoadField(TString s){ |
void TrkTrack::LoadField(TString path){ |
290 |
readb_(s.Data()); |
|
291 |
|
strcpy(path_.path,path.Data()); |
292 |
|
path_.pathlen = path.Length(); |
293 |
|
path_.error = 0; |
294 |
|
readb_(); |
295 |
|
|
296 |
}; |
}; |
297 |
/** |
/** |
298 |
* Tracking method. It calls F77 mini routine. |
* Tracking method. It calls F77 mini routine. |
303 |
|
|
304 |
extern cMini2track track_; |
extern cMini2track track_; |
305 |
fail = 0; |
fail = 0; |
|
// extern cMini2fitinfo fit_info_; |
|
|
// extern void mini_2_(int*,int*); |
|
|
// track_.xm[0]=1.0; |
|
|
|
|
|
cout << "fitting routine "<<endl; |
|
306 |
|
|
307 |
for(int i=0; i<6; i++) track_.xm[i]=xm[i]; |
for(int i=0; i<6; i++) track_.xm[i]=xm[i]; |
308 |
for(int i=0; i<6; i++) track_.ym[i]=ym[i]; |
for(int i=0; i<6; i++) track_.ym[i]=ym[i]; |
314 |
|
|
315 |
// initial guess of "al" with linear fit |
// initial guess of "al" with linear fit |
316 |
// if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.){ |
// if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.){ |
|
|
|
317 |
// cout << "initial guess "<<endl; |
// cout << "initial guess "<<endl; |
|
|
|
318 |
// double szz=0., szx=0., szy=0., ssx=0., ssy=0., sz=0., s1=0.; |
// double szz=0., szx=0., szy=0., ssx=0., ssy=0., sz=0., s1=0.; |
319 |
// double det, ax, ay, bx, by; |
// double det, ax, ay, bx, by; |
320 |
// for(int i=0; i<NPLANE; i++) { |
// for(int i=0; i<NPLANE; i++) { |
340 |
// if(ax<0.) al[3]=acos(-1.)-al[3]; |
// if(ax<0.) al[3]=acos(-1.)-al[3]; |
341 |
// } |
// } |
342 |
// al[4]=0.; |
// al[4]=0.; |
|
|
|
343 |
// } |
// } |
344 |
// end guess |
// end guess |
345 |
|
|
346 |
for(int i=0; i<5; i++) track_.al[i]=al[i]; |
for(int i=0; i<5; i++) track_.al[i]=al[i]; |
347 |
track_.zini = 23.5; // ZINI = 23.5 !!! it should be the same parameter in all codes |
track_.zini = 23.5; |
348 |
|
// ZINI = 23.5 !!! it should be the same parameter in all codes |
349 |
|
|
|
// elena ----- |
|
350 |
if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.)guess_(); |
if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.)guess_(); |
|
// elena ----- |
|
351 |
|
|
352 |
|
// --------------------- free momentum |
353 |
if(pfixed==0.) { |
if(pfixed==0.) { |
354 |
// al[4]=0.; // free momentum |
// al[4]=0.; // free momentum |
355 |
track_.pfixed=0.; // " |
track_.pfixed=0.; // " |
356 |
} |
} |
357 |
|
// --------------------- fixed momentum |
358 |
if(pfixed!=0.) { |
if(pfixed!=0.) { |
359 |
al[4]=1./pfixed; // to fix the momentum |
al[4]=1./pfixed; // to fix the momentum |
360 |
track_.pfixed=pfixed; // " |
track_.pfixed=pfixed; // " |
979 |
* Loads the magnetic field. |
* Loads the magnetic field. |
980 |
* @param s Path of the magnetic-field files. |
* @param s Path of the magnetic-field files. |
981 |
*/ |
*/ |
982 |
void TrkLevel2::LoadField(TString s){ |
void TrkLevel2::LoadField(TString path){ |
983 |
readb_(s.Data()); |
// |
984 |
|
strcpy(path_.path,path.Data()); |
985 |
|
path_.pathlen = path.Length(); |
986 |
|
path_.error = 0; |
987 |
|
readb_(); |
988 |
|
// |
989 |
}; |
}; |
990 |
//-------------------------------------- |
//-------------------------------------- |
991 |
// |
// |