1 /*
2 * PropZZ_Test.java
3 *
4 * Created on July 24, 2007, 11:02 PM
5 *
6 * $Id: PropZZ_Test.java,v 1.1.1.1 2010/04/08 20:38:00 jeremy Exp $
7 */
8
9 package org.lcsim.recon.tracking.trfzp;
10
11 import junit.framework.TestCase;
12 import org.lcsim.recon.tracking.spacegeom.SpacePath;
13 import org.lcsim.recon.tracking.spacegeom.SpacePoint;
14 import org.lcsim.recon.tracking.trfbase.ETrack;
15 import org.lcsim.recon.tracking.trfbase.PropDir;
16 import org.lcsim.recon.tracking.trfbase.PropStat;
17 import org.lcsim.recon.tracking.trfbase.Propagator;
18 import org.lcsim.recon.tracking.trfbase.Surface;
19 import org.lcsim.recon.tracking.trfbase.TrackDerivative;
20 import org.lcsim.recon.tracking.trfbase.TrackError;
21 import org.lcsim.recon.tracking.trfbase.TrackVector;
22 import org.lcsim.recon.tracking.trfbase.VTrack;
23 import org.lcsim.recon.tracking.trfutil.Assert;
24
25 /**
26 *
27 * @author Norman Graf
28 */
29 public class PropZZ_Test extends TestCase
30 {
31 private boolean debug;
32 /** Creates a new instance of PropZZ_Test */
33 public void testPropZZ()
34 {
35 String ok_prefix = "PropZZ (I): ";
36 String error_prefix = "PropZZ test (E): ";
37
38 if(debug) System.out.println( ok_prefix
39 + "-------- Testing component PropZZ. --------" );
40
41 {
42 if(debug) System.out.println( ok_prefix + "Test constructor." );
43 PropZZ prop = new PropZZ(2.0);
44 if(debug) System.out.println( prop );
45
46 //********************************************************************
47
48 // Here we propagate some tracks both forward and backward and then
49 // each back to the original track. We check that the returned
50 // track parameters match those of the original track.
51 if(debug) System.out.println( ok_prefix + "Check reversibility." );
52
53 double z1[] ={10., 10., 10., 10., 10., 10.};
54 double z2[] ={20., -20., 20., 20., -20., -20.};
55 double z3[] ={30., -30., 30., 30., -30., -30.};
56 int sign_dz[]={ 1, -1, 1, 1, -1, -1 };
57 double x[] ={ -2., 2., 40., 40., -2., -2. };
58 double y[] ={ 3., 3., 3., 3., 3., 3. };
59 double dxdz[]={-1.5, -2.3, 0., 1.5, 0., 0. };
60 double dydz[]={ 2.3, -1.5, -2.3, 0., 2.3, 0. };
61 double qp[] ={ 0.05, -0.05, 0.05, -0.05, -0.05, 0.05 };
62
63 double maxdiff = 1.e-12;
64 int ntrk = 6;
65 int i;
66 for ( i=0; i<ntrk; ++i )
67 {
68 if(debug) System.out.println( "********** Propagate track " + i + ". **********" );
69 PropStat pstat = new PropStat();
70 SurfZPlane szp1 = new SurfZPlane(z1[i]);
71 SurfZPlane szp2 = new SurfZPlane(z2[i]);
72 SurfZPlane szp3 = new SurfZPlane(z3[i]);
73 TrackVector vec1 = new TrackVector();
74 vec1.set(SurfZPlane.IX , x[i]); // x
75 vec1.set(SurfZPlane.IY , y[i]); // y
76 vec1.set(SurfZPlane.IDXDZ, dxdz[i]); // dx/dz
77 vec1.set(SurfZPlane.IDYDZ, dydz[i]); // dy/dz
78 vec1.set(SurfZPlane.IQP , qp[i]); // q/p
79
80 VTrack trv1 = new VTrack(szp2.newPureSurface(),vec1);
81 if (sign_dz[i]==1)
82 {trv1.setForward();
83 }
84 else
85 {trv1.setBackward();
86 }
87 if(debug) System.out.println( " starting: " + trv1 );
88 VTrack trv2f = new VTrack(trv1);
89 pstat = prop.vecDirProp(trv2f,szp3,PropDir.FORWARD);
90 Assert.assertTrue( pstat.forward() );
91 if(debug) System.out.println( " forward: " + trv2f );
92 VTrack trv2b = new VTrack(trv1);
93 pstat = prop.vecDirProp(trv2b,szp1,PropDir.BACKWARD);
94 Assert.assertTrue( pstat.backward() );
95 if(debug) System.out.println( " backward: " + trv2b );
96 VTrack trv2fb = new VTrack(trv2f);
97 pstat = prop.vecDirProp(trv2fb,szp2,PropDir.BACKWARD);
98 Assert.assertTrue( pstat.backward() );
99 if(debug) System.out.println( " f return: " + trv2fb );
100 VTrack trv2bf = new VTrack(trv2b);
101 pstat = prop.vecDirProp(trv2bf,szp2,PropDir.FORWARD);
102 Assert.assertTrue( pstat.forward() );
103 if(debug) System.out.println( " b return: " + trv2bf );
104 double difff =
105 szp1.vecDiff(trv2fb.vector(),trv1.vector()).amax();
106 double diffb =
107 szp1.vecDiff(trv2bf.vector(),trv1.vector()).amax();
108 if(debug) System.out.println( "diffs: " + difff + ' ' + diffb );
109 Assert.assertTrue( difff < maxdiff );
110 Assert.assertTrue( diffb < maxdiff );
111 }
112
113 // Repeat the above with errors.
114 if(debug) System.out.println( ok_prefix + "Check reversibility with errors." );
115 double exx[] = { 0.01, 0.01, 0.01, 0.01, 0.01, 0.01 };
116 double exy[] = { 0.01, -0.01, 0.01, -0.01, 0.01, -0.01 };
117 double eyy[] = { 0.25, 0.25, 0.25, 0.25, 0.25, 0.25, };
118 double exdx[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
119 double eydx[] = { 0.04, -0.04, 0.04, -0.04, 0.04, -0.04, };
120 double edxdx[] = { 0.01, 0.01, 0.01, 0.01, 0.01, 0.01 };
121 double exdy[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
122 double edxdy[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
123 double eydy[] = { 0.04, -0.04, 0.04, -0.04, 0.04, -0.04 };
124 double edydy[] = { 0.02, 0.02, 0.02, 0.02, 0.02, 0.02 };
125 double exqp[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
126 double eyqp[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
127 double edxqp[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
128 double edyqp[] = { 0.004, -0.004, 0.004, -0.004, 0.004, -0.004 };
129 double eqpqp[] = { 0.01, 0.01, 0.01, 0.01, 0.01, 0.01 };
130
131 for ( i=0; i<ntrk; ++i )
132 {
133 if(debug) System.out.println( "********** Propagate track " + i + ". **********" );
134 PropStat pstat = new PropStat();
135 SurfZPlane szp1= new SurfZPlane((z1[i]));
136 SurfZPlane szp2= new SurfZPlane((z2[i]));
137 SurfZPlane szp3= new SurfZPlane((z3[i]));
138 TrackVector vec1 = new TrackVector();
139 vec1.set(SurfZPlane.IX , x[i]); // x
140 vec1.set(SurfZPlane.IY , y[i]); // y
141 vec1.set(SurfZPlane.IDXDZ , dxdz[i]); // dx/dz
142 vec1.set(SurfZPlane.IDYDZ , dydz[i]); // dy/dz
143 vec1.set(SurfZPlane.IQP , qp[i]); // q/p
144
145 TrackError err1 = new TrackError();
146 err1.set(SurfZPlane.IX,SurfZPlane.IX , exx[i]);
147 err1.set(SurfZPlane.IX,SurfZPlane.IY , exy[i]);
148 err1.set(SurfZPlane.IY,SurfZPlane.IY , eyy[i]);
149 err1.set(SurfZPlane.IX,SurfZPlane.IDXDZ , exdx[i]);
150 err1.set(SurfZPlane.IY,SurfZPlane.IDXDZ , eydx[i]);
151 err1.set(SurfZPlane.IDXDZ,SurfZPlane.IDXDZ , edxdx[i]);
152 err1.set(SurfZPlane.IX,SurfZPlane.IDYDZ , exdy[i]);
153 err1.set(SurfZPlane.IY,SurfZPlane.IDYDZ , eydy[i]);
154 err1.set(SurfZPlane.IDXDZ,SurfZPlane.IDYDZ , edxdy[i]);
155 err1.set(SurfZPlane.IDYDZ,SurfZPlane.IDYDZ , edydy[i]);
156 err1.set(SurfZPlane.IX,SurfZPlane.IQP , exqp[i]);
157 err1.set(SurfZPlane.IY,SurfZPlane.IQP , eyqp[i]);
158 err1.set(SurfZPlane.IDXDZ,SurfZPlane.IQP , edxqp[i]);
159 err1.set(SurfZPlane.IDYDZ,SurfZPlane.IQP , edyqp[i]);
160 err1.set(SurfZPlane.IQP,SurfZPlane.IQP , eqpqp[i]);
161 ETrack trv1 = new ETrack(szp2.newPureSurface(),vec1,err1);
162 if (sign_dz[i]==1)
163 {trv1.setForward();
164 }
165 else
166 {trv1.setBackward();
167 }
168 if(debug) System.out.println( " starting: " + trv1 );
169
170 ETrack trv2f = new ETrack(trv1);
171 pstat = prop.errDirProp(trv2f,szp3,PropDir.FORWARD);
172 Assert.assertTrue( pstat.forward() );
173 if(debug) System.out.println( " forward: " + trv2f );
174 ETrack trv2b = new ETrack(trv1);
175 pstat = prop.errDirProp(trv2b,szp1,PropDir.BACKWARD);
176 Assert.assertTrue( pstat.backward() );
177 if(debug) System.out.println( " backward: " + trv2b );
178 ETrack trv2fb = new ETrack(trv2f);
179 pstat = prop.errDirProp(trv2fb,szp2,PropDir.BACKWARD);
180 Assert.assertTrue( pstat.backward() );
181 if(debug) System.out.println( " f return: " + trv2fb );
182 ETrack trv2bf = new ETrack(trv2b);
183 pstat = prop.errDirProp(trv2bf,szp2,PropDir.FORWARD);
184 Assert.assertTrue( pstat.forward() );
185 if(debug) System.out.println( " b return: " + trv2bf );
186 double difff =
187 szp1.vecDiff(trv2fb.vector(),trv1.vector()).amax();
188 double diffb =
189 szp1.vecDiff(trv2bf.vector(),trv1.vector()).amax();
190 if(debug) System.out.println( "vec diffs: " + difff + ' ' + diffb );
191 Assert.assertTrue( difff < maxdiff );
192 Assert.assertTrue( diffb < maxdiff );
193
194 TrackError dfb = trv2fb.error().minus(trv1.error());
195 TrackError dbf = trv2bf.error().minus(trv1.error());
196 double edifff = dfb.amax();
197 double ediffb = dbf.amax();
198 if(debug) System.out.println( "err diffs: " + edifff + ' ' + ediffb );
199 Assert.assertTrue( edifff < maxdiff );
200 Assert.assertTrue( ediffb < maxdiff );
201
202 }
203
204 //********************************************************************
205
206 if(debug) System.out.println( ok_prefix + "Test Nearest Propagation" );
207
208 PropStat pstat = new PropStat();
209
210 SurfZPlane szp1 = new SurfZPlane(5.);
211 SurfZPlane szp2 = new SurfZPlane(10.);
212
213 TrackVector vec1 = new TrackVector();
214 vec1.set(SurfZPlane.IX , 1.); // x
215 vec1.set(SurfZPlane.IY , 1.); // y
216 vec1.set(SurfZPlane.IDXDZ , 1.); // dx/dz
217 vec1.set(SurfZPlane.IDYDZ , 1.); // dy/dz
218 vec1.set(SurfZPlane.IQP , 0.01); // q/p
219
220 VTrack trv1 = new VTrack(szp1.newPureSurface(),vec1);
221 trv1.setForward();
222
223 if(debug) System.out.println( " starting: " + trv1 );
224 VTrack trv2n = new VTrack(trv1);
225 pstat = prop.vecDirProp(trv2n,szp2,PropDir.NEAREST_MOVE);
226 Assert.assertTrue( pstat.forward() );
227 if(debug) System.out.println( " nearest: " + trv2n );
228
229 trv1.setBackward();
230
231 if(debug) System.out.println( " starting: " + trv1 );
232 trv2n = new VTrack(trv1);
233 pstat = prop.vecDirProp(trv2n,szp2,PropDir.NEAREST);
234 Assert.assertTrue( pstat.backward() );
235 if(debug) System.out.println( " nearest: " + trv2n );
236
237 //********************************************************************
238
239 if(debug) System.out.println( ok_prefix + "Test same surface" );
240 VTrack trvs0 = new VTrack(szp1.newPureSurface(),vec1);
241 trvs0.setForward();
242 VTrack trvs = new VTrack(trvs0);
243 PropStat tst = prop.vecDirProp(trvs,szp1,PropDir.NEAREST_MOVE);
244 Assert.assertTrue( !tst.success() );
245 tst = prop.vecDirProp(trvs,szp1,PropDir.FORWARD_MOVE);
246 Assert.assertTrue( !tst.success() );
247 tst = prop.vecDirProp(trvs,szp1,PropDir.BACKWARD_MOVE);
248 Assert.assertTrue( !tst.success() );
249 tst = prop.vecDirProp(trvs,szp1,PropDir.NEAREST);
250 Assert.assertTrue( tst.success() );
251 tst = prop.vecDirProp(trvs,szp1,PropDir.FORWARD);
252 Assert.assertTrue( tst.success() );
253 tst = prop.vecDirProp(trvs,szp1,PropDir.BACKWARD);
254 Assert.assertTrue( tst.success() );
255 Assert.assertTrue( trvs0.equals(trvs) );
256
257 //********************************************************************
258
259 if(debug) System.out.println( ok_prefix + "Test cloning." );
260 Assert.assertTrue( prop.newPropagator() != null);
261
262
263 //********************************************************************
264
265 if(debug) System.out.println( ok_prefix + "Test the field." );
266 Assert.assertTrue( prop.bField() == 2.0 );
267
268 //********************************************************************
269 }
270
271 if(debug) System.out.println( ok_prefix + "Test Zero B field propagation" );
272
273 {
274 PropZZ prop0 = new PropZZ(0.);
275 if(debug) System.out.println( prop0 );
276 Assert.assertTrue( prop0.bField() == 0. );
277
278 double z=10.;
279 Surface srf = new SurfZPlane(z);
280 VTrack trv0 = new VTrack(srf);
281 TrackVector vec = new TrackVector();
282 vec.set(SurfZPlane.IX, 2.);
283 vec.set(SurfZPlane.IY, 10.);
284 vec.set(SurfZPlane.IDXDZ, 4.);
285 vec.set(SurfZPlane.IDYDZ, 2.);
286 trv0.setVector(vec);
287 trv0.setForward();
288 z=4.;
289 Surface srf_to = new SurfZPlane(z);
290
291 VTrack trv = new VTrack(trv0);
292 PropStat pstat = prop0.vecDirProp(trv,srf_to,PropDir.FORWARD);
293 Assert.assertTrue( !pstat.success() );
294
295 trv = new VTrack(trv0);
296 pstat = prop0.vecDirProp(trv,srf_to,PropDir.BACKWARD);
297 Assert.assertTrue( pstat.success() );
298
299 trv = new VTrack(trv0);
300 trv.setBackward();
301 pstat = prop0.vecDirProp(trv,srf_to,PropDir.BACKWARD);
302 Assert.assertTrue( !pstat.success() );
303
304 trv = new VTrack(trv0);
305 trv.setBackward();
306 pstat = prop0.vecDirProp(trv,srf_to,PropDir.FORWARD);
307 Assert.assertTrue( pstat.success() );
308
309 trv = new VTrack(trv0);
310 trv.setForward();
311 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
312 Assert.assertTrue( pstat.success() );
313
314 Assert.assertTrue( pstat.backward() );
315 Assert.assertTrue( trv.vector(SurfZPlane.IDXDZ) == trv0.vector(SurfZPlane.IDXDZ) );
316 Assert.assertTrue( trv.vector(SurfZPlane.IDYDZ) == trv0.vector(SurfZPlane.IDYDZ) );
317 Assert.assertTrue(trv.surface().pureEqual(srf_to));
318
319 check_zero_propagation(trv0,trv,pstat);
320
321 srf_to = new SurfZPlane(6.);
322 trv = new VTrack(trv0);
323 trv.setForward();
324 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
325 Assert.assertTrue( pstat.success() );
326 check_zero_propagation(trv0,trv,pstat);
327
328 srf_to = new SurfZPlane(14.);
329 trv = new VTrack(trv0);
330 trv.setForward();
331 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
332 Assert.assertTrue( pstat.success() );
333 check_zero_propagation(trv0,trv,pstat);
334
335 srf_to = new SurfZPlane(-1.);
336 trv = new VTrack(trv0);
337 trv.setForward();
338 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
339 Assert.assertTrue( pstat.success() );
340 check_zero_propagation(trv0,trv,pstat);
341
342 srf_to = new SurfZPlane(-14.);
343 trv = new VTrack(trv0);
344 trv.setForward();
345 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
346 Assert.assertTrue( pstat.success() );
347 check_zero_propagation(trv0,trv,pstat);
348
349 srf_to = new SurfZPlane(14.);
350 trv = new VTrack(trv0);
351 trv.setSurface( new SurfZPlane(1.));
352 trv.setBackward();
353 VTrack tmp = new VTrack(trv);
354 VTrack der = new VTrack(trv);
355 VTrack der1 = new VTrack(trv);
356 der1.setForward();
357 pstat = prop0.vecDirProp(trv,srf_to,PropDir.NEAREST);
358 Assert.assertTrue( pstat.success() );
359 check_zero_propagation(tmp,trv,pstat);
360 check_derivatives(prop0,der,srf_to);
361 check_derivatives(prop0,der1,srf_to);
362
363 }
364
365
366 if(debug) System.out.println( ok_prefix
367 + "------------- All tests passed. -------------" );
368
369 //********************************************************************
370 }
371
372
373 private static void check_zero_propagation( VTrack trv0, VTrack trv, PropStat pstat)
374 {
375
376 SpacePoint sp = trv.spacePoint();
377 SpacePoint sp0 = trv0.spacePoint();
378
379 SpacePath sv = trv.spacePath();
380 SpacePath sv0 = trv0.spacePath();
381
382 Assert.assertTrue( Math.abs(sv0.dx() - sv.dx())<1e-7 );
383 Assert.assertTrue( Math.abs(sv0.dy() - sv.dy())<1e-7 );
384 Assert.assertTrue( Math.abs(sv0.dz() - sv.dz())<1e-7 );
385
386 double x0 = sp0.x();
387 double y0 = sp0.y();
388 double z0 = sp0.z();
389 double x1 = sp.x();
390 double y1 = sp.y();
391 double z1 = sp.z();
392
393 double dx = sv.dx();
394 double dy = sv.dy();
395 double dz = sv.dz();
396
397 double prod = dx*(x1-x0)+dy*(y1-y0)+dz*(z1-z0);
398 double moda = Math.sqrt((x1-x0)*(x1-x0)+(y1-y0)*(y1-y0) + (z1-z0)*(z1-z0));
399 double modb = Math.sqrt(dx*dx+dy*dy+dz*dz);
400 double st = pstat.pathDistance();
401 Assert.assertTrue( Math.abs(prod-st) < 1.e-7 );
402 Assert.assertTrue( Math.abs(Math.abs(prod) - moda*modb) < 1.e-7 );
403 }
404
405 private static void check_derivatives( Propagator prop, VTrack trv0, Surface srf)
406 {
407 for(int i=0;i<4;++i)
408 for(int j=0;j<4;++j)
409 check_derivative(prop,trv0,srf,i,j);
410 }
411
412 private static void check_derivative( Propagator prop, VTrack trv0, Surface srf,int i,int j)
413 {
414
415 double dx = 1.e-3;
416 VTrack trv = new VTrack(trv0);
417 TrackVector vec = trv.vector();
418 boolean forward = trv0.isForward();
419
420 VTrack trv_0 = new VTrack(trv0);
421 TrackDerivative der = new TrackDerivative();
422 PropStat pstat = prop.vecProp(trv_0,srf,der);
423 Assert.assertTrue(pstat.success());
424
425 TrackVector tmp = new TrackVector(vec);
426 tmp.set(j, tmp.get(j)+dx);
427 trv.setVector(tmp);
428 if(forward) trv.setForward();
429 else trv.setBackward();
430
431 VTrack trv_pl = new VTrack(trv);
432 pstat = prop.vecProp(trv_pl,srf);
433 Assert.assertTrue(pstat.success());
434
435 TrackVector vecpl = trv_pl.vector();
436
437 tmp= new TrackVector(vec);
438 tmp.set(j, tmp.get(j)-dx);
439 trv.setVector(tmp);
440 if(forward) trv.setForward();
441 else trv.setBackward();
442
443 VTrack trv_mn = new VTrack(trv);
444 pstat = prop.vecProp(trv_mn,srf);
445 Assert.assertTrue(pstat.success());
446
447 TrackVector vecmn = trv_mn.vector();
448
449 double dy = (vecpl.get(i)-vecmn.get(i))/2.;
450
451 double dydx = dy/dx;
452
453 double didj = der.get(i,j);
454
455 if( Math.abs(didj) > 1e-10 )
456 Assert.assertTrue( Math.abs((dydx - didj)/didj) < 1e-4 );
457 else
458 Assert.assertTrue( Math.abs(dydx) < 1e-4 );
459 }
460 }