| 18 |
include '../common/level1.f' |
include '../common/level1.f' |
| 19 |
include '../common/level2.f' |
include '../common/level2.f' |
| 20 |
|
|
| 21 |
include '../common/momanhough_init.f' |
c include '../common/momanhough_init.f' |
| 22 |
|
|
| 23 |
* flag to set debug mode |
* flag to set debug mode |
| 24 |
logical DEBUG |
logical DEBUG |
| 61 |
# endif |
# endif |
| 62 |
|
|
| 63 |
integer minevent !first event to be analysed |
integer minevent !first event to be analysed |
| 64 |
logical FIMAGE ! |
c logical FIMAGE ! |
| 65 |
|
|
| 66 |
|
|
| 67 |
COMMON/QUEST/IQUEST(100) |
COMMON/QUEST/IQUEST(100) |
| 406 |
print*,iev,' ** ',nev2 |
print*,iev,' ** ',nev2 |
| 407 |
endif |
endif |
| 408 |
|
|
| 409 |
*------------------------------------------------------------------------------- |
* /////////////////////////////////////////////// |
| 410 |
* STEP 1 |
* \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\ |
|
*------------------------------------------------------------------------------- |
|
|
* X-Y cluster association |
|
|
* |
|
|
* Clusters are associated to form COUPLES |
|
|
* Clusters not associated in any couple are called SINGLETS |
|
|
* |
|
|
* Track identification (Hough transform) and fitting is first done on couples. |
|
|
* Hence singlets are possibly added to the track. |
|
|
* |
|
|
* Variables assigned by the routine "cl_to_couples" are those in the |
|
|
* common blocks: |
|
|
* - common/clusters/cl_good |
|
|
* - common/couples/clx,cly,ncp_plane,ncp_tot,cp_useds1,cp_useds2 |
|
|
* - common/singlets/ncls,cls,cl_single |
|
|
*------------------------------------------------------------------------------- |
|
|
*------------------------------------------------------------------------------- |
|
| 411 |
|
|
| 412 |
iflag=0 |
iflag=0 |
| 413 |
call cl_to_couples(iflag) |
call track_finding(iflag) |
| 414 |
if(iflag.eq.1)then !bad event |
if(iflag.eq.1)then !bad event |
| 415 |
goto 880 !fill ntp and go to next event |
goto 880 !fill ntp and go to next event |
| 416 |
endif |
endif |
| 417 |
|
c$$$*------------------------------------------------------------------------------- |
| 418 |
*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
c$$$* STEP 1 |
| 419 |
* selezione di tracce pulite per diagnostica |
c$$$*------------------------------------------------------------------------------- |
| 420 |
*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
c$$$* X-Y cluster association |
| 421 |
c$$$ if(DEBUG)then |
c$$$* |
| 422 |
c$$$ do ip=1,nplanes |
c$$$* Clusters are associated to form COUPLES |
| 423 |
c$$$ if(ncp_plane(ip).ne.1)good2=.false. |
c$$$* Clusters not associated in any couple are called SINGLETS |
| 424 |
c$$$ enddo |
c$$$* |
| 425 |
c$$$c if(good2.eq.0)goto 100!next event |
c$$$* Track identification (Hough transform) and fitting is first done on couples. |
| 426 |
c$$$c if(good2.eq.0)goto 880!fill ntp and go to next event |
c$$$* Hence singlets are possibly added to the track. |
| 427 |
|
c$$$* |
| 428 |
|
c$$$* Variables assigned by the routine "cl_to_couples" are those in the |
| 429 |
|
c$$$* common blocks: |
| 430 |
|
c$$$* - common/clusters/cl_good |
| 431 |
|
c$$$* - common/couples/clx,cly,ncp_plane,ncp_tot,cp_useds1,cp_useds2 |
| 432 |
|
c$$$* - common/singlets/ncls,cls,cl_single |
| 433 |
|
c$$$*------------------------------------------------------------------------------- |
| 434 |
|
c$$$*------------------------------------------------------------------------------- |
| 435 |
|
c$$$ |
| 436 |
|
c$$$ iflag=0 |
| 437 |
|
c$$$ call cl_to_couples(iflag) |
| 438 |
|
c$$$ if(iflag.eq.1)then !bad event |
| 439 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 440 |
c$$$ endif |
c$$$ endif |
| 441 |
|
c$$$ |
| 442 |
|
c$$$*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 443 |
|
c$$$* selezione di tracce pulite per diagnostica |
| 444 |
*----------------------------------------------------- |
c$$$*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 445 |
*----------------------------------------------------- |
c$$$c$$$ if(DEBUG)then |
| 446 |
* HOUGH TRASFORM |
c$$$c$$$ do ip=1,nplanes |
| 447 |
*----------------------------------------------------- |
c$$$c$$$ if(ncp_plane(ip).ne.1)good2=.false. |
| 448 |
*----------------------------------------------------- |
c$$$c$$$ enddo |
| 449 |
|
c$$$c$$$c if(good2.eq.0)goto 100!next event |
| 450 |
|
c$$$c$$$c if(good2.eq.0)goto 880!fill ntp and go to next event |
| 451 |
*------------------------------------------------------------------------------- |
c$$$c$$$ endif |
| 452 |
* STEP 2 |
c$$$ |
| 453 |
*------------------------------------------------------------------------------- |
c$$$ |
| 454 |
* |
c$$$ |
| 455 |
* Association of couples to form |
c$$$*----------------------------------------------------- |
| 456 |
* - DOUBLETS in YZ view |
c$$$*----------------------------------------------------- |
| 457 |
* - TRIPLETS in XZ view |
c$$$* HOUGH TRASFORM |
| 458 |
* |
c$$$*----------------------------------------------------- |
| 459 |
* Variables assigned by the routine "cp_to_doubtrip" are those in the |
c$$$*----------------------------------------------------- |
| 460 |
* common blocks: |
c$$$ |
| 461 |
* - common/hough_param/ |
c$$$ |
| 462 |
* $ alfayz1, !Y0 |
c$$$*------------------------------------------------------------------------------- |
| 463 |
* $ alfayz2, !tg theta-yz |
c$$$* STEP 2 |
| 464 |
* $ alfaxz1, !X0 |
c$$$*------------------------------------------------------------------------------- |
| 465 |
* $ alfaxz2, !tg theta-xz |
c$$$* |
| 466 |
* $ alfaxz3 !1/r |
c$$$* Association of couples to form |
| 467 |
* - common/doublets/ndblt,cpyz1,cpyz2 |
c$$$* - DOUBLETS in YZ view |
| 468 |
* - common/triplets/ntrpt,cpxz1,cpxz2,cpxz3 |
c$$$* - TRIPLETS in XZ view |
| 469 |
*------------------------------------------------------------------------------- |
c$$$* |
| 470 |
*------------------------------------------------------------------------------- |
c$$$* Variables assigned by the routine "cp_to_doubtrip" are those in the |
| 471 |
|
c$$$* common blocks: |
| 472 |
iflag=0 |
c$$$* - common/hough_param/ |
| 473 |
call cp_to_doubtrip(iflag) |
c$$$* $ alfayz1, !Y0 |
| 474 |
if(iflag.eq.1)then !bad event |
c$$$* $ alfayz2, !tg theta-yz |
| 475 |
goto 880 !fill ntp and go to next event |
c$$$* $ alfaxz1, !X0 |
| 476 |
endif |
c$$$* $ alfaxz2, !tg theta-xz |
| 477 |
|
c$$$* $ alfaxz3 !1/r |
| 478 |
|
c$$$* - common/doublets/ndblt,cpyz1,cpyz2 |
| 479 |
*------------------------------------------------------------------------------- |
c$$$* - common/triplets/ntrpt,cpxz1,cpxz2,cpxz3 |
| 480 |
* STEP 3 |
c$$$*------------------------------------------------------------------------------- |
| 481 |
*------------------------------------------------------------------------------- |
c$$$*------------------------------------------------------------------------------- |
| 482 |
* |
c$$$ |
| 483 |
* Classification of doublets and triplets to form CLOUDS, |
c$$$ iflag=0 |
| 484 |
* according to distance in parameter space. |
c$$$ call cp_to_doubtrip(iflag) |
| 485 |
* |
c$$$ if(iflag.eq.1)then !bad event |
| 486 |
* cloud = cluster of points (doublets/triplets) in parameter space |
c$$$ goto 880 !fill ntp and go to next event |
| 487 |
* |
c$$$ endif |
| 488 |
* |
c$$$ |
| 489 |
* |
c$$$ |
| 490 |
* Variables assigned by the routine "doub_to_YZcloud" are those in the |
c$$$*------------------------------------------------------------------------------- |
| 491 |
* common blocks: |
c$$$* STEP 3 |
| 492 |
* - common/clouds_yz/ |
c$$$*------------------------------------------------------------------------------- |
| 493 |
* $ nclouds_yz |
c$$$* |
| 494 |
* $ ,alfayz1_av,alfayz2_av |
c$$$* Classification of doublets and triplets to form CLOUDS, |
| 495 |
* $ ,ptcloud_yz,db_cloud,cpcloud_yz |
c$$$* according to distance in parameter space. |
| 496 |
* |
c$$$* |
| 497 |
* Variables assigned by the routine "trip_to_XZcloud" are those in the |
c$$$* cloud = cluster of points (doublets/triplets) in parameter space |
| 498 |
* common blocks: |
c$$$* |
| 499 |
* common/clouds_xz/ |
c$$$* |
| 500 |
* $ nclouds_xz xz2_av,alfaxz3_av |
c$$$* |
| 501 |
* $ ,ptcloud_xz,tr_cloud,cpcloud_xz |
c$$$* Variables assigned by the routine "doub_to_YZcloud" are those in the |
| 502 |
*------------------------------------------------------------------------------- |
c$$$* common blocks: |
| 503 |
*------------------------------------------------------------------------------- |
c$$$* - common/clouds_yz/ |
| 504 |
|
c$$$* $ nclouds_yz |
| 505 |
iflag=0 |
c$$$* $ ,alfayz1_av,alfayz2_av |
| 506 |
call doub_to_YZcloud(iflag) |
c$$$* $ ,ptcloud_yz,db_cloud,cpcloud_yz |
| 507 |
if(iflag.eq.1)then !bad event |
c$$$* |
| 508 |
goto 880 !fill ntp and go to next event |
c$$$* Variables assigned by the routine "trip_to_XZcloud" are those in the |
| 509 |
endif |
c$$$* common blocks: |
| 510 |
iflag=0 |
c$$$* common/clouds_xz/ |
| 511 |
call trip_to_XZcloud(iflag) |
c$$$* $ nclouds_xz xz2_av,alfaxz3_av |
| 512 |
if(iflag.eq.1)then !bad event |
c$$$* $ ,ptcloud_xz,tr_cloud,cpcloud_xz |
| 513 |
goto 880 !fill ntp and go to next event |
c$$$*------------------------------------------------------------------------------- |
| 514 |
endif |
c$$$*------------------------------------------------------------------------------- |
| 515 |
|
c$$$ |
| 516 |
|
c$$$ iflag=0 |
| 517 |
|
c$$$ call doub_to_YZcloud(iflag) |
| 518 |
|
c$$$ if(iflag.eq.1)then !bad event |
| 519 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 520 |
|
c$$$ endif |
| 521 |
|
c$$$ iflag=0 |
| 522 |
|
c$$$ call trip_to_XZcloud(iflag) |
| 523 |
|
c$$$ if(iflag.eq.1)then !bad event |
| 524 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 525 |
|
c$$$ endif |
| 526 |
|
c$$$ |
| 527 |
c***************************************************** |
c***************************************************** |
| 528 |
cccccc 01/12/2005 modified by elena |
cccccc 01/12/2005 modified by elena |
| 529 |
if(DEBUG)then |
if(DEBUG)then |
| 534 |
endif |
endif |
| 535 |
c***************************************************** |
c***************************************************** |
| 536 |
|
|
|
*------------------------------------------------------------------------------- |
|
|
* STEP 4 (ITERATED until any other physical track isn't found) |
|
|
*------------------------------------------------------------------------------- |
|
|
* |
|
|
* YZ and XZ clouds are combined in order to obtain the initial guess |
|
|
* of the candidate-track parameters. |
|
|
* A minimum number of matching couples between YZ and XZ clouds is required. |
|
|
* |
|
|
* A TRACK CANDIDATE is defined by |
|
|
* - the couples resulting from the INTERSECTION of the two clouds, and |
|
|
* - the associated track parameters (evaluated by performing a zero-order |
|
|
* track fitting) |
|
|
* |
|
|
* The NTRACKS candidate-track parameters are stored in common block: |
|
|
* |
|
|
* - common/track_candidates/NTRACKS,AL_STORE |
|
|
* $ ,XV_STORE,YV_STORE,ZV_STORE |
|
|
* $ ,XM_STORE,YM_STORE,ZM_STORE |
|
|
* $ ,RESX_STORE,RESY_STORE |
|
|
* $ ,AXV_STORE,AYV_STORE |
|
|
* $ ,XGOOD_STORE,YGOOD_STORE |
|
|
* $ ,CP_STORE,RCHI2_STORE |
|
|
* |
|
|
*------------------------------------------------------------------------------- |
|
|
*------------------------------------------------------------------------------- |
|
|
ntrk=0 !counter of identified physical tracks |
|
|
|
|
|
11111 continue !<<<<<<< come here when performing a new search |
|
| 537 |
|
|
| 538 |
iflag=0 |
iflag=0 |
| 539 |
call clouds_to_ctrack(iflag) |
call track_fitting(iflag) |
|
if(iflag.eq.1)then !no candidate tracks found |
|
|
goto 880 !fill ntp and go to next event |
|
|
endif |
|
|
|
|
|
FIMAGE=.false. !processing best track (not track image) |
|
|
ibest=0 !best track among candidates |
|
|
iimage=0 !track image |
|
|
* ------------- select the best track ------------- |
|
|
rchi2best=1000000000. |
|
|
do i=1,ntracks |
|
|
if(RCHI2_STORE(i).lt.rchi2best.and. |
|
|
$ RCHI2_STORE(i).gt.0)then |
|
|
ibest=i |
|
|
rchi2best=RCHI2_STORE(i) |
|
|
endif |
|
|
enddo |
|
|
if(ibest.eq.0)goto 880 !>> no good candidates |
|
|
*------------------------------------------------------------------------------- |
|
|
* The best track candidate (ibest) is selected and a new fitting is performed. |
|
|
* Previous to this, the track is refined by: |
|
|
* - possibly adding new COUPLES or SINGLETS from the missing planes |
|
|
* - evaluating the coordinates with improved PFAs |
|
|
* ( angle-dependent ETA algorithms ) |
|
|
*------------------------------------------------------------------------------- |
|
|
|
|
|
1212 continue !<<<<< come here to fit track-image |
|
|
|
|
|
if(.not.FIMAGE)then !processing best candidate |
|
|
icand=ibest |
|
|
else !processing image |
|
|
icand=iimage |
|
|
iimage=0 |
|
|
endif |
|
|
if(icand.eq.0)then |
|
|
print*,'HAI FATTO UN CASINO!!!!!! icand = ',icand |
|
|
$ ,ibest,iimage |
|
|
return |
|
|
endif |
|
|
|
|
|
* *-*-*-*-*-*-*-*-*-*-*-*-*-*-* |
|
|
call refine_track(icand) |
|
|
* *-*-*-*-*-*-*-*-*-*-*-*-*-*-* |
|
|
|
|
|
* ********************************************************** |
|
|
* ************************** FIT *** FIT *** FIT *** FIT *** |
|
|
* ********************************************************** |
|
|
do i=1,5 |
|
|
AL(i)=dble(AL_STORE(i,icand)) |
|
|
enddo |
|
|
ifail=0 !error flag in chi2 computation |
|
|
jstep=0 !# minimization steps |
|
|
|
|
|
call mini_2(jstep,ifail) |
|
|
if(ifail.ne.0) then |
|
|
if(DEBUG)then |
|
|
print *, |
|
|
$ '*** MINIMIZATION FAILURE *** (mini_2) ' |
|
|
$ ,iev |
|
|
endif |
|
|
chi2=-chi2 |
|
|
endif |
|
|
|
|
|
if(DEBUG)then |
|
|
print*,'----------------------------- improved track coord' |
|
|
22222 format(i2,' * ',3f10.4,' --- ',4f10.4,' --- ',2f4.0,2f10.5) |
|
|
do ip=1,6 |
|
|
write(*,22222)ip,zm(ip),xm(ip),ym(ip) |
|
|
$ ,xm_A(ip),ym_A(ip),xm_B(ip),ym_B(ip) |
|
|
$ ,xgood(ip),ygood(ip),resx(ip),resy(ip) |
|
|
enddo |
|
|
endif |
|
|
|
|
|
c rchi2=chi2/dble(ndof) |
|
|
if(DEBUG)then |
|
|
print*,' ' |
|
|
print*,'****** SELECTED TRACK *************' |
|
|
print*,'# R. chi2 RIG' |
|
|
print*,' --- ',chi2,' --- ' |
|
|
$ ,1./abs(AL(5)) |
|
|
print*,'***********************************' |
|
|
endif |
|
|
* ********************************************************** |
|
|
* ************************** FIT *** FIT *** FIT *** FIT *** |
|
|
* ********************************************************** |
|
|
|
|
|
|
|
|
* ------------- search if the track has an IMAGE ------------- |
|
|
* ------------- (also this is stored ) ------------- |
|
|
if(FIMAGE)goto 122 !>>> jump! (this is already an image) |
|
|
* now search for track-image, by comparing couples IDs |
|
|
do i=1,ntracks |
|
|
iimage=i |
|
|
do ip=1,nplanes |
|
|
if( CP_STORE(nplanes-ip+1,icand).ne. |
|
|
$ -1*CP_STORE(nplanes-ip+1,i) )iimage=0 |
|
|
enddo |
|
|
if( iimage.ne.0.and. |
|
|
c $ RCHI2_STORE(i).le.CHI2MAX.and. |
|
|
c $ RCHI2_STORE(i).gt.0.and. |
|
|
$ .true.)then |
|
|
if(DEBUG)print*,'Track candidate ',iimage |
|
|
$ ,' >>> TRACK IMAGE >>> of' |
|
|
$ ,ibest |
|
|
goto 122 !image track found |
|
|
endif |
|
|
enddo |
|
|
122 continue |
|
|
|
|
|
* --- and store the results -------------------------------- |
|
|
ntrk = ntrk + 1 !counter of found tracks |
|
|
if(.not.FIMAGE |
|
|
$ .and.iimage.eq.0) image(ntrk)= 0 |
|
|
if(.not.FIMAGE |
|
|
$ .and.iimage.ne.0)image(ntrk)=ntrk+1 !this is the image of the next |
|
|
if(FIMAGE) image(ntrk)=ntrk-1 !this is the image of the previous |
|
|
|
|
|
call fill_level2_tracks(ntrk) !==> good2=.true. |
|
|
c print*,'++++++++++ iimage,fimage,ntrk,image ' |
|
|
c $ ,iimage,fimage,ntrk,image(ntrk) |
|
|
|
|
|
if(ntrk.eq.NTRKMAX)then |
|
|
if(DEBUG) |
|
|
$ print*, |
|
|
$ '** warning ** number of identified '// |
|
|
$ 'tracks exceeds vector dimension ' |
|
|
$ ,'( ',NTRKMAX,' )' |
|
|
cc good2=.false. |
|
|
goto 880 !fill ntp and go to next event |
|
|
endif |
|
|
if(iimage.ne.0)then |
|
|
FIMAGE=.true. ! |
|
|
goto 1212 !>>> fit image-track |
|
|
endif |
|
|
|
|
|
* --- then remove selected clusters (ibest+iimage) from clouds ---- |
|
|
call clean_XYclouds(ibest) |
|
| 540 |
if(iflag.eq.1)then !bad event |
if(iflag.eq.1)then !bad event |
| 541 |
goto 880 !fill ntp and go to next event |
goto 880 !fill ntp and go to next event |
| 542 |
endif |
endif |
| 543 |
|
|
|
* ********************************************************** |
|
|
* condition to start a new search |
|
|
* ********************************************************** |
|
|
ixznew=0 |
|
|
do ixz=1,nclouds_xz |
|
|
if(ptcloud_xz(ixz).ge.nptxz_min)ixznew=1 |
|
|
enddo |
|
|
iyznew=0 |
|
|
do iyz=1,nclouds_yz |
|
|
if(ptcloud_yz(iyz).ge.nptyz_min)iyznew=1 |
|
|
enddo |
|
|
|
|
|
if(ixznew.ne.0.and. |
|
|
$ iyznew.ne.0.and. |
|
|
$ rchi2best.le.CHI2MAX.and. |
|
|
c $ rchi2best.lt.15..and. |
|
|
$ .true.)then |
|
|
if(DEBUG)then |
|
|
print*,'***** NEW SEARCH ****' |
|
|
endif |
|
|
goto 11111 !try new search |
|
|
|
|
|
endif |
|
|
* ********************************************** |
|
| 544 |
|
|
| 545 |
|
c$$$ |
| 546 |
|
c$$$*------------------------------------------------------------------------------- |
| 547 |
|
c$$$* STEP 4 (ITERATED until any other physical track isn't found) |
| 548 |
|
c$$$*------------------------------------------------------------------------------- |
| 549 |
|
c$$$* |
| 550 |
|
c$$$* YZ and XZ clouds are combined in order to obtain the initial guess |
| 551 |
|
c$$$* of the candidate-track parameters. |
| 552 |
|
c$$$* A minimum number of matching couples between YZ and XZ clouds is required. |
| 553 |
|
c$$$* |
| 554 |
|
c$$$* A TRACK CANDIDATE is defined by |
| 555 |
|
c$$$* - the couples resulting from the INTERSECTION of the two clouds, and |
| 556 |
|
c$$$* - the associated track parameters (evaluated by performing a zero-order |
| 557 |
|
c$$$* track fitting) |
| 558 |
|
c$$$* |
| 559 |
|
c$$$* The NTRACKS candidate-track parameters are stored in common block: |
| 560 |
|
c$$$* |
| 561 |
|
c$$$* - common/track_candidates/NTRACKS,AL_STORE |
| 562 |
|
c$$$* $ ,XV_STORE,YV_STORE,ZV_STORE |
| 563 |
|
c$$$* $ ,XM_STORE,YM_STORE,ZM_STORE |
| 564 |
|
c$$$* $ ,RESX_STORE,RESY_STORE |
| 565 |
|
c$$$* $ ,AXV_STORE,AYV_STORE |
| 566 |
|
c$$$* $ ,XGOOD_STORE,YGOOD_STORE |
| 567 |
|
c$$$* $ ,CP_STORE,RCHI2_STORE |
| 568 |
|
c$$$* |
| 569 |
|
c$$$*------------------------------------------------------------------------------- |
| 570 |
|
c$$$*------------------------------------------------------------------------------- |
| 571 |
|
c$$$ ntrk=0 !counter of identified physical tracks |
| 572 |
|
c$$$ |
| 573 |
|
c$$$11111 continue !<<<<<<< come here when performing a new search |
| 574 |
|
c$$$ |
| 575 |
|
c$$$ iflag=0 |
| 576 |
|
c$$$ call clouds_to_ctrack(iflag) |
| 577 |
|
c$$$ if(iflag.eq.1)then !no candidate tracks found |
| 578 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 579 |
|
c$$$ endif |
| 580 |
|
c$$$ |
| 581 |
|
c$$$ FIMAGE=.false. !processing best track (not track image) |
| 582 |
|
c$$$ ibest=0 !best track among candidates |
| 583 |
|
c$$$ iimage=0 !track image |
| 584 |
|
c$$$* ------------- select the best track ------------- |
| 585 |
|
c$$$ rchi2best=1000000000. |
| 586 |
|
c$$$ do i=1,ntracks |
| 587 |
|
c$$$ if(RCHI2_STORE(i).lt.rchi2best.and. |
| 588 |
|
c$$$ $ RCHI2_STORE(i).gt.0)then |
| 589 |
|
c$$$ ibest=i |
| 590 |
|
c$$$ rchi2best=RCHI2_STORE(i) |
| 591 |
|
c$$$ endif |
| 592 |
|
c$$$ enddo |
| 593 |
|
c$$$ if(ibest.eq.0)goto 880 !>> no good candidates |
| 594 |
|
c$$$*------------------------------------------------------------------------------- |
| 595 |
|
c$$$* The best track candidate (ibest) is selected and a new fitting is performed. |
| 596 |
|
c$$$* Previous to this, the track is refined by: |
| 597 |
|
c$$$* - possibly adding new COUPLES or SINGLETS from the missing planes |
| 598 |
|
c$$$* - evaluating the coordinates with improved PFAs |
| 599 |
|
c$$$* ( angle-dependent ETA algorithms ) |
| 600 |
|
c$$$*------------------------------------------------------------------------------- |
| 601 |
|
c$$$ |
| 602 |
|
c$$$ 1212 continue !<<<<< come here to fit track-image |
| 603 |
|
c$$$ |
| 604 |
|
c$$$ if(.not.FIMAGE)then !processing best candidate |
| 605 |
|
c$$$ icand=ibest |
| 606 |
|
c$$$ else !processing image |
| 607 |
|
c$$$ icand=iimage |
| 608 |
|
c$$$ iimage=0 |
| 609 |
|
c$$$ endif |
| 610 |
|
c$$$ if(icand.eq.0)then |
| 611 |
|
c$$$ print*,'HAI FATTO UN CASINO!!!!!! icand = ',icand |
| 612 |
|
c$$$ $ ,ibest,iimage |
| 613 |
|
c$$$ return |
| 614 |
|
c$$$ endif |
| 615 |
|
c$$$ |
| 616 |
|
c$$$* *-*-*-*-*-*-*-*-*-*-*-*-*-*-* |
| 617 |
|
c$$$ call refine_track(icand) |
| 618 |
|
c$$$* *-*-*-*-*-*-*-*-*-*-*-*-*-*-* |
| 619 |
|
c$$$ |
| 620 |
|
c$$$* ********************************************************** |
| 621 |
|
c$$$* ************************** FIT *** FIT *** FIT *** FIT *** |
| 622 |
|
c$$$* ********************************************************** |
| 623 |
|
c$$$ do i=1,5 |
| 624 |
|
c$$$ AL(i)=dble(AL_STORE(i,icand)) |
| 625 |
|
c$$$ enddo |
| 626 |
|
c$$$ ifail=0 !error flag in chi2 computation |
| 627 |
|
c$$$ jstep=0 !# minimization steps |
| 628 |
|
c$$$ |
| 629 |
|
c$$$ call mini_2(jstep,ifail) |
| 630 |
|
c$$$ if(ifail.ne.0) then |
| 631 |
|
c$$$ if(DEBUG)then |
| 632 |
|
c$$$ print *, |
| 633 |
|
c$$$ $ '*** MINIMIZATION FAILURE *** (mini_2) ' |
| 634 |
|
c$$$ $ ,iev |
| 635 |
|
c$$$ endif |
| 636 |
|
c$$$ chi2=-chi2 |
| 637 |
|
c$$$ endif |
| 638 |
|
c$$$ |
| 639 |
|
c$$$ if(DEBUG)then |
| 640 |
|
c$$$ print*,'----------------------------- improved track coord' |
| 641 |
|
c$$$22222 format(i2,' * ',3f10.4,' --- ',4f10.4,' --- ',2f4.0,2f10.5) |
| 642 |
|
c$$$ do ip=1,6 |
| 643 |
|
c$$$ write(*,22222)ip,zm(ip),xm(ip),ym(ip) |
| 644 |
|
c$$$ $ ,xm_A(ip),ym_A(ip),xm_B(ip),ym_B(ip) |
| 645 |
|
c$$$ $ ,xgood(ip),ygood(ip),resx(ip),resy(ip) |
| 646 |
|
c$$$ enddo |
| 647 |
|
c$$$ endif |
| 648 |
|
c$$$ |
| 649 |
|
c$$$c rchi2=chi2/dble(ndof) |
| 650 |
|
c$$$ if(DEBUG)then |
| 651 |
|
c$$$ print*,' ' |
| 652 |
|
c$$$ print*,'****** SELECTED TRACK *************' |
| 653 |
|
c$$$ print*,'# R. chi2 RIG' |
| 654 |
|
c$$$ print*,' --- ',chi2,' --- ' |
| 655 |
|
c$$$ $ ,1./abs(AL(5)) |
| 656 |
|
c$$$ print*,'***********************************' |
| 657 |
|
c$$$ endif |
| 658 |
|
c$$$* ********************************************************** |
| 659 |
|
c$$$* ************************** FIT *** FIT *** FIT *** FIT *** |
| 660 |
|
c$$$* ********************************************************** |
| 661 |
|
c$$$ |
| 662 |
|
c$$$ |
| 663 |
|
c$$$* ------------- search if the track has an IMAGE ------------- |
| 664 |
|
c$$$* ------------- (also this is stored ) ------------- |
| 665 |
|
c$$$ if(FIMAGE)goto 122 !>>> jump! (this is already an image) |
| 666 |
|
c$$$* now search for track-image, by comparing couples IDs |
| 667 |
|
c$$$ do i=1,ntracks |
| 668 |
|
c$$$ iimage=i |
| 669 |
|
c$$$ do ip=1,nplanes |
| 670 |
|
c$$$ if( CP_STORE(nplanes-ip+1,icand).ne. |
| 671 |
|
c$$$ $ -1*CP_STORE(nplanes-ip+1,i) )iimage=0 |
| 672 |
|
c$$$ enddo |
| 673 |
|
c$$$ if( iimage.ne.0.and. |
| 674 |
|
c$$$c $ RCHI2_STORE(i).le.CHI2MAX.and. |
| 675 |
|
c$$$c $ RCHI2_STORE(i).gt.0.and. |
| 676 |
|
c$$$ $ .true.)then |
| 677 |
|
c$$$ if(DEBUG)print*,'Track candidate ',iimage |
| 678 |
|
c$$$ $ ,' >>> TRACK IMAGE >>> of' |
| 679 |
|
c$$$ $ ,ibest |
| 680 |
|
c$$$ goto 122 !image track found |
| 681 |
|
c$$$ endif |
| 682 |
|
c$$$ enddo |
| 683 |
|
c$$$ 122 continue |
| 684 |
|
c$$$ |
| 685 |
|
c$$$* --- and store the results -------------------------------- |
| 686 |
|
c$$$ ntrk = ntrk + 1 !counter of found tracks |
| 687 |
|
c$$$ if(.not.FIMAGE |
| 688 |
|
c$$$ $ .and.iimage.eq.0) image(ntrk)= 0 |
| 689 |
|
c$$$ if(.not.FIMAGE |
| 690 |
|
c$$$ $ .and.iimage.ne.0)image(ntrk)=ntrk+1 !this is the image of the next |
| 691 |
|
c$$$ if(FIMAGE) image(ntrk)=ntrk-1 !this is the image of the previous |
| 692 |
|
c$$$ |
| 693 |
|
c$$$ call fill_level2_tracks(ntrk) !==> good2=.true. |
| 694 |
|
c$$$c print*,'++++++++++ iimage,fimage,ntrk,image ' |
| 695 |
|
c$$$c $ ,iimage,fimage,ntrk,image(ntrk) |
| 696 |
|
c$$$ |
| 697 |
|
c$$$ if(ntrk.eq.NTRKMAX)then |
| 698 |
|
c$$$ if(DEBUG) |
| 699 |
|
c$$$ $ print*, |
| 700 |
|
c$$$ $ '** warning ** number of identified '// |
| 701 |
|
c$$$ $ 'tracks exceeds vector dimension ' |
| 702 |
|
c$$$ $ ,'( ',NTRKMAX,' )' |
| 703 |
|
c$$$cc good2=.false. |
| 704 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 705 |
|
c$$$ endif |
| 706 |
|
c$$$ if(iimage.ne.0)then |
| 707 |
|
c$$$ FIMAGE=.true. ! |
| 708 |
|
c$$$ goto 1212 !>>> fit image-track |
| 709 |
|
c$$$ endif |
| 710 |
|
c$$$ |
| 711 |
|
c$$$* --- then remove selected clusters (ibest+iimage) from clouds ---- |
| 712 |
|
c$$$ call clean_XYclouds(ibest) |
| 713 |
|
c$$$ if(iflag.eq.1)then !bad event |
| 714 |
|
c$$$ goto 880 !fill ntp and go to next event |
| 715 |
|
c$$$ endif |
| 716 |
|
c$$$ |
| 717 |
|
c$$$* ********************************************************** |
| 718 |
|
c$$$* condition to start a new search |
| 719 |
|
c$$$* ********************************************************** |
| 720 |
|
c$$$ ixznew=0 |
| 721 |
|
c$$$ do ixz=1,nclouds_xz |
| 722 |
|
c$$$ if(ptcloud_xz(ixz).ge.nptxz_min)ixznew=1 |
| 723 |
|
c$$$ enddo |
| 724 |
|
c$$$ iyznew=0 |
| 725 |
|
c$$$ do iyz=1,nclouds_yz |
| 726 |
|
c$$$ if(ptcloud_yz(iyz).ge.nptyz_min)iyznew=1 |
| 727 |
|
c$$$ enddo |
| 728 |
|
c$$$ |
| 729 |
|
c$$$ if(ixznew.ne.0.and. |
| 730 |
|
c$$$ $ iyznew.ne.0.and. |
| 731 |
|
c$$$ $ rchi2best.le.CHI2MAX.and. |
| 732 |
|
c$$$c $ rchi2best.lt.15..and. |
| 733 |
|
c$$$ $ .true.)then |
| 734 |
|
c$$$ if(DEBUG)then |
| 735 |
|
c$$$ print*,'***** NEW SEARCH ****' |
| 736 |
|
c$$$ endif |
| 737 |
|
c$$$ goto 11111 !try new search |
| 738 |
|
c$$$ |
| 739 |
|
c$$$ endif |
| 740 |
|
c$$$* ********************************************** |
| 741 |
|
c$$$ |
| 742 |
* >>>>>>>>>>>>>>>>>>> NT-UPLE filling <<<<<<<<<<<<<<<<<<<< |
* >>>>>>>>>>>>>>>>>>> NT-UPLE filling <<<<<<<<<<<<<<<<<<<< |
| 743 |
* + + + + + + + + + + + + + + + + + |
* + + + + + + + + + + + + + + + + + |
| 744 |
* / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / |
* / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / \ / |
| 789 |
|
|
| 790 |
8800 continue |
8800 continue |
| 791 |
|
|
| 792 |
|
|
| 793 |
|
* /////////////////////////////////////////////// |
| 794 |
|
* \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\ |
| 795 |
|
|
| 796 |
call HCDIR('//LEVEL2',' ') |
call HCDIR('//LEVEL2',' ') |
| 797 |
call HFNT(ntp_level2) !fill LEVEL2 nt-uple |
call HFNT(ntp_level2) !fill LEVEL2 nt-uple |
| 798 |
if(.not.DBG_FILLED.and.DEBUG) |
if(.not.DBG_FILLED.and.DEBUG) |
Parent Directory
| 
Revision Log
| 
Patch