| 1 |
************************************************************************* |
************************************************************************ |
| 2 |
* |
* |
| 3 |
* functions.f |
* functions.f |
| 4 |
* |
* |
| 65 |
|
|
| 66 |
include 'commontracker.f' |
include 'commontracker.f' |
| 67 |
|
|
| 68 |
|
view = view |
| 69 |
nld=INT((istrip-1)/nstrips_ladder)+1 |
nld=INT((istrip-1)/nstrips_ladder)+1 |
| 70 |
|
|
| 71 |
end |
end |
| 132 |
|
|
| 133 |
end |
end |
| 134 |
|
|
|
c$$$ |
|
|
c$$$c$$$c------------------------------------------------------------------------ |
|
|
c$$$c$$$ |
|
|
c$$$c$$$c NB: le coordinate in mech_pos.dat sono calcolate a partire da alcuni dati |
|
|
c$$$c$$$c contenuti in commontracker.f. forse si puo' evitare mech_pos.dat e mettere |
|
|
c$$$c$$$c tutto in commontracker.f |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ subroutine mech_sensor !it reads sensors coordinates (in PAMELA reference |
|
|
c$$$c$$$ ! frame) from a text file and it uses them to fill |
|
|
c$$$c$$$ ! x/y/z_mech_sensor variables, taking into account |
|
|
c$$$c$$$ ! last plane inversion |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ include './commontracker.f' |
|
|
c$$$c$$$ include './common_tracks.f' |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ real xvec(nladders_view),yvec(2),zvec(nplanes) |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ integer id !file identifier |
|
|
c$$$c$$$ logical od !.true. if the specified unit is connected to a file |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ do id=20,100,1 !opens the file using a free file id |
|
|
c$$$c$$$ inquire (id, opened=od) |
|
|
c$$$c$$$ if(.not.od) goto 666 |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$ 666 continue |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ open(id,FILE='../common/mech_pos.dat') !sensors centres coordinates in mm in |
|
|
c$$$c$$$ ! PAMELA reference frame: |
|
|
c$$$c$$$ ! the first plane is the one with lowest Z (the one |
|
|
c$$$c$$$ ! nearest the calorimeter) |
|
|
c$$$c$$$ ! the first ladder is the one with lowest X (the |
|
|
c$$$c$$$ ! one on which the first X strip is) |
|
|
c$$$c$$$ ! the first sensor is the one with lowest Y (the |
|
|
c$$$c$$$ ! one on which the first Y strip is) for planes |
|
|
c$$$c$$$ ! 2..6. for plane 1 the first sensor has higher Y |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ read(20,*) xvec |
|
|
c$$$c$$$ read(20,*) yvec |
|
|
c$$$c$$$ read(20,*) zvec |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ do i=1,nplanes |
|
|
c$$$c$$$ do j=1,nladders_view |
|
|
c$$$c$$$ do k=1,2 |
|
|
c$$$c$$$ x_mech_sensor(i,j,k)=xvec(j) |
|
|
c$$$c$$$ y_mech_sensor(i,j,k)=yvec(k) |
|
|
c$$$c$$$ z_mech_sensor(i,j,k)=zvec(i) |
|
|
c$$$c$$$ if(i.eq.1) then !y coordinates of first plane (11th view) are |
|
|
c$$$c$$$ y_mech_sensor(i,j,k)=-yvec(k) ! exchanged due to last plane inversion |
|
|
c$$$c$$$ endif |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ close(id) |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ |
|
|
c$$$c$$$c$$$ ! *** INIZIO DEBUG *** |
|
|
c$$$c$$$c$$$ do i=1,6 |
|
|
c$$$c$$$c$$$ do j=1,3 |
|
|
c$$$c$$$c$$$ do k=1,2 |
|
|
c$$$c$$$c$$$ c print*,x_mech_sensor(1,j,k) |
|
|
c$$$c$$$c$$$ print*,y_mech_sensor(i,j,k) |
|
|
c$$$c$$$c$$$ c print*,z_mech_sensor(i,j,k) |
|
|
c$$$c$$$c$$$ enddo |
|
|
c$$$c$$$c$$$ enddo |
|
|
c$$$c$$$c$$$ print*,' ' |
|
|
c$$$c$$$c$$$ enddo |
|
|
c$$$c$$$c$$$ ! *** FINE DEBUG *** |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ |
|
|
c$$$c$$$ return |
|
|
c$$$c$$$ end |
|
|
|
|
|
|
|
|
c$$$c------------------------------------------------------------------------ |
|
|
c$$$ |
|
|
c$$$c NB: le coordinate in mech_pos.dat sono calcolate a partire da alcuni dati |
|
|
c$$$c contenuti in commontracker.f. forse si puo' evitare mech_pos.dat e mettere |
|
|
c$$$c tutto in commontracker.f |
|
|
c$$$ |
|
|
c$$$ |
|
|
c$$$ subroutine mech_sensor |
|
|
c$$$c !it reads sensors coordinates (in PAMELA reference |
|
|
c$$$c ! frame) from a text file and it uses them to fill |
|
|
c$$$c ! x/y/z_mech_sensor variables, taking into account |
|
|
c$$$c ! last plane inversion |
|
|
c$$$ |
|
|
c$$$ include './commontracker.f' |
|
|
c$$$ include './common_tracks.f' |
|
|
c$$$ |
|
|
c$$$ real xvec(nladders_view),yvec(2),zvec(nplanes) |
|
|
c$$$ |
|
|
c$$$ integer id !file identifier |
|
|
c$$$ logical od !.true. if the specified unit is connected to a file |
|
|
c$$$ |
|
|
c$$$ do id=20,100,1 !opens the file using a free file id |
|
|
c$$$ inquire (id, opened=od) |
|
|
c$$$ if(.not.od) goto 666 |
|
|
c$$$ enddo |
|
|
c$$$ 666 continue |
|
|
c$$$ |
|
|
c$$$c open(id,FILE='../common/mech_pos.dat') !sensors centres coordinates in mm in |
|
|
c$$$c open(id,FILE='source/common/mech_pos.dat') |
|
|
c$$$c call system('cp $TRK_GRND/source/common/mech_pos.dat .') |
|
|
c$$$ print *,'Opening file: mech_pos.dat' |
|
|
c$$$ open(id,FILE='./bin-aux/mech_pos.dat',IOSTAT=iostat) |
|
|
c$$$c !sensors centres coordinates in mm in |
|
|
c$$$c ! PAMELA reference frame: |
|
|
c$$$c ! the first plane is the one with lowest Z (the one |
|
|
c$$$c ! nearest the calorimeter) |
|
|
c$$$c ! the first ladder is the one with lowest X (the |
|
|
c$$$c ! one on which the first X strip is) |
|
|
c$$$c ! the first sensor is the one with lowest Y (the |
|
|
c$$$c ! one on which the first Y strip is) for planes |
|
|
c$$$c ! 2..6. for plane 1 the first sensor has higher Y |
|
|
c$$$ |
|
|
c$$$ if(iostat.ne.0)then |
|
|
c$$$ print*,'MECH_SENSOR: *** Error in opening file ***' |
|
|
c$$$ return |
|
|
c$$$ endif |
|
|
c$$$ |
|
|
c$$$ read(id,*) xvec |
|
|
c$$$ read(id,*) yvec |
|
|
c$$$ read(id,*) zvec |
|
|
c$$$ |
|
|
c$$$ do i=1,nplanes |
|
|
c$$$ do j=1,nladders_view |
|
|
c$$$ do k=1,2 |
|
|
c$$$ x_mech_sensor(i,j,k)=xvec(j) |
|
|
c$$$ y_mech_sensor(i,j,k)=yvec(k) |
|
|
c$$$ z_mech_sensor(i,j,k)=zvec(i) |
|
|
c$$$ if(i.eq.1) then !y coordinates of first plane (11th view) are |
|
|
c$$$ y_mech_sensor(i,j,k)=-yvec(k) ! exchanged due to last plane inversion |
|
|
c$$$ endif |
|
|
c$$$ enddo |
|
|
c$$$ enddo |
|
|
c$$$ enddo |
|
|
c$$$ |
|
|
c$$$ close(id) |
|
|
c$$$c call system('rm -f mech_pos.dat') |
|
|
c$$$ |
|
|
c$$$c$$$ ! *** INIZIO DEBUG *** |
|
|
c$$$c$$$ do i=1,6 |
|
|
c$$$c$$$c do j=1,3 |
|
|
c$$$c$$$ do k=1,2 |
|
|
c$$$c$$$ j=1 |
|
|
c$$$c$$$c print*,x_mech_sensor(1,j,k) |
|
|
c$$$c$$$ print*,y_mech_sensor(i,j,k) |
|
|
c$$$c$$$c print*,z_mech_sensor(i,j,k) |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$c enddo |
|
|
c$$$c$$$ print*,' ' |
|
|
c$$$c$$$ enddo |
|
|
c$$$c$$$ ! *** FINE DEBUG *** |
|
|
c$$$ |
|
|
c$$$ |
|
|
c$$$ return |
|
|
c$$$ end |
|
|
c$$$ |
|
|
c$$$ |
|
|
c$$$c------------------------------------------------------------------------ |
|
| 135 |
|
|
| 136 |
|
c------------------------------------------------------------------------ |
| 137 |
|
|
| 138 |
function coordsi(istrip,view) |
function coordsi(istrip,view) |
| 139 |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
| 164 |
|
|
| 165 |
if(mod(view,2).eq.0) then !X view |
if(mod(view,2).eq.0) then !X view |
| 166 |
|
|
| 167 |
if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
| 168 |
print*,'functions: WARNING: false X strip: strip ',is |
c print*,'functions: WARNING: false X strip: strip ',is |
| 169 |
endif |
c endif |
| 170 |
|
|
| 171 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
| 172 |
|
|
| 235 |
|
|
| 236 |
if(mod(view,2).eq.0) then !X view |
if(mod(view,2).eq.0) then !X view |
| 237 |
|
|
| 238 |
if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
| 239 |
print*,'functions: WARNING: false X strip: strip ',is |
c print*,'functions: WARNING: false X strip: strip ',is |
| 240 |
endif |
c endif |
| 241 |
|
|
| 242 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
| 243 |
|
|
| 275 |
c------------------------------------------------------------------------ |
c------------------------------------------------------------------------ |
| 276 |
|
|
| 277 |
|
|
| 278 |
|
double precision function dcoordsi(rstrip,view) |
| 279 |
|
* |
| 280 |
|
* same as COORDSI, but accept a real value of strip!!! |
| 281 |
|
* and gives a double precision output |
| 282 |
|
* |
| 283 |
|
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
| 284 |
|
* it gives the strip coordinate in micrometers, |
| 285 |
|
* knowing the strip number (1..3072) and the view |
| 286 |
|
* number. the origin of the coordinate is on the |
| 287 |
|
* centre of the sensor the strip belongs to. |
| 288 |
|
* the axes directions are the same as in the PAMELA |
| 289 |
|
* reference frame (i.e.: the 11th view coordinate |
| 290 |
|
* direction has to be inverted here) |
| 291 |
|
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
| 292 |
|
|
| 293 |
|
c integer is,view,istrip |
| 294 |
|
|
| 295 |
|
integer view,is,istrip |
| 296 |
|
real rstrip |
| 297 |
|
double precision strip,stripladder,p |
| 298 |
|
double precision edge,dim |
| 299 |
|
double precision coord1 |
| 300 |
|
|
| 301 |
|
|
| 302 |
|
include 'commontracker.f' |
| 303 |
|
|
| 304 |
|
c NB mettere il 1024 nel commontracker...!??? |
| 305 |
|
|
| 306 |
|
strip=DBLE(rstrip) |
| 307 |
|
|
| 308 |
|
istrip = int(strip+0.5) !istrip stores the closest integer to strip |
| 309 |
|
|
| 310 |
|
is=istrip !it stores istrip number |
| 311 |
|
is=mod(is-1,1024)+1 !it puts all clusters on a single ladder |
| 312 |
|
|
| 313 |
|
dcoordsi=0. |
| 314 |
|
|
| 315 |
|
if(mod(view,2).eq.0) then !X view |
| 316 |
|
|
| 317 |
|
c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
| 318 |
|
c print*,'functions: WARNING: false X strip: strip ',is |
| 319 |
|
c endif |
| 320 |
|
|
| 321 |
|
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
| 322 |
|
|
| 323 |
|
edge=edgeX |
| 324 |
|
dim=SiDimX |
| 325 |
|
|
| 326 |
|
elseif(mod(view,2).eq.1) then !Y view |
| 327 |
|
|
| 328 |
|
edge=edgeY |
| 329 |
|
dim=SiDimY |
| 330 |
|
|
| 331 |
|
c$$$ if(view.eq.11) then !INVERSIONE!??? |
| 332 |
|
c$$$ is=1025-is |
| 333 |
|
c$$$ endif |
| 334 |
|
|
| 335 |
|
endif |
| 336 |
|
|
| 337 |
|
|
| 338 |
|
stripladder = DBLE(is)+(strip-DBLE(istrip))!cluster position relative to ladder |
| 339 |
|
p=pitch(view) |
| 340 |
|
|
| 341 |
|
ccccc coord1=(is-1)*p !referred to 1st sensor strip |
| 342 |
|
coord1=(stripladder-1)*p !referred to 1st sensor strip |
| 343 |
|
coord1=coord1+edge !referred to sensor edge |
| 344 |
|
dcoordsi=coord1-dim/2 !referred to the centre of the sensor |
| 345 |
|
|
| 346 |
|
if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views |
| 347 |
|
dcoordsi=-dcoordsi |
| 348 |
|
endif |
| 349 |
|
|
| 350 |
|
end |
| 351 |
|
|
| 352 |
|
|
| 353 |
|
|
| 354 |
|
c------------------------------------------------------------------------ |
| 355 |
|
|
| 356 |
|
|
| 357 |
function coord(coordsi,view,ladder,sen) |
function coord(coordsi,view,ladder,sen) |
| 358 |
* it gives the coordinate in |
* it gives the coordinate in |
| 359 |
* micrometers, knowing the coordinate in the sensor |
* micrometers, knowing the coordinate in the sensor |
| 438 |
|
|
| 439 |
dcoord=0. |
dcoord=0. |
| 440 |
|
|
| 441 |
|
if( |
| 442 |
|
$ ladder.lt.1.or. |
| 443 |
|
$ ladder.gt.nladders_view.or. |
| 444 |
|
$ sen.lt.1.or. |
| 445 |
|
$ sen.gt.2.or. |
| 446 |
|
$ view.lt.1.or. |
| 447 |
|
$ view.gt.nviews.or. |
| 448 |
|
$ .false.)then |
| 449 |
|
print*,'dcoord ---> wrong input: ladder ',ladder |
| 450 |
|
$ ,' sensor ',sen |
| 451 |
|
$ ,' view ',view |
| 452 |
|
return |
| 453 |
|
endif |
| 454 |
|
|
| 455 |
sx=ladder |
sx=ladder |
| 456 |
sy=sen |
sy=sen |
| 457 |
sz=npl(view) |
sz=npl(view) |
| 458 |
|
|
| 459 |
|
|
| 460 |
if(mod(view,2).eq.0) then !X view |
if(mod(view,2).eq.0) then !X view |
| 461 |
|
|
| 462 |
trasl=x_mech_sensor(sz,sx,sy) !in mm |
trasl=x_mech_sensor(sz,sx,sy) !in mm |
| 477 |
|
|
| 478 |
|
|
| 479 |
c------------------------------------------------------------------------ |
c------------------------------------------------------------------------ |
| 480 |
|
integer function nsatstrips(ic) |
| 481 |
|
*-------------------------------------------------------------- |
| 482 |
|
* this function returns the number of saturated strips |
| 483 |
|
* inside a cluster |
| 484 |
|
*-------------------------------------------------------------- |
| 485 |
|
include 'commontracker.f' |
| 486 |
|
include 'level1.f' |
| 487 |
|
include 'calib.f' |
| 488 |
|
|
| 489 |
|
|
| 490 |
|
|
| 491 |
|
integer nsat |
| 492 |
|
nsat = 0 |
| 493 |
|
iv=VIEW(ic) |
| 494 |
|
if(mod(iv,2).eq.1)incut=incuty |
| 495 |
|
if(mod(iv,2).eq.0)incut=incutx |
| 496 |
|
istart = INDSTART(IC) |
| 497 |
|
istop = TOTCLLENGTH |
| 498 |
|
if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1 |
| 499 |
|
do i = INDMAX(IC),istart,-1 |
| 500 |
|
c cut = incut*CLSIGMA(i) |
| 501 |
|
c if(CLSIGNAL(i).ge.cut)then |
| 502 |
|
if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx) |
| 503 |
|
$ .or. |
| 504 |
|
$ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then |
| 505 |
|
nsat = nsat +1 |
| 506 |
|
else |
| 507 |
|
goto 10 |
| 508 |
|
endif |
| 509 |
|
enddo |
| 510 |
|
10 continue |
| 511 |
|
do i = INDMAX(IC)+1,istop |
| 512 |
|
c cut = incut*CLSIGMA(i) |
| 513 |
|
c if(CLSIGNAL(i).ge.cut)then |
| 514 |
|
if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx) |
| 515 |
|
$ .or. |
| 516 |
|
$ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then |
| 517 |
|
nsat = nsat +1 |
| 518 |
|
else |
| 519 |
|
goto 20 |
| 520 |
|
endif |
| 521 |
|
enddo |
| 522 |
|
20 continue |
| 523 |
|
|
| 524 |
|
nsatstrips = nsat |
| 525 |
|
return |
| 526 |
|
end |
| 527 |
|
|
| 528 |
|
c------------------------------------------------------------------------ |
| 529 |
|
integer function nbadstrips(ncog,ic) |
| 530 |
|
*-------------------------------------------------------------- |
| 531 |
|
* this function returns the number of BAD strips |
| 532 |
|
* inside a cluster: |
| 533 |
|
* - if NCOG=0, the number BAD strips inside the whole cluster |
| 534 |
|
* are given, according to the cluster multiplicity |
| 535 |
|
* |
| 536 |
|
* - if NCOG>0, the number BAD strips is evaluated using NCOG |
| 537 |
|
* strips, even if they have a negative signal (according to Landi) |
| 538 |
|
*-------------------------------------------------------------- |
| 539 |
|
include 'commontracker.f' |
| 540 |
|
include 'level1.f' |
| 541 |
|
include 'calib.f' |
| 542 |
|
|
| 543 |
|
integer nbad |
| 544 |
|
nbad = 0 |
| 545 |
|
|
| 546 |
|
if (ncog.gt.0) then |
| 547 |
|
|
| 548 |
|
* --> signal of the central strip |
| 549 |
|
sc = CLSIGNAL(INDMAX(ic)) !center |
| 550 |
|
* signal of adjacent strips |
| 551 |
|
sl1 = -100000 !left 1 |
| 552 |
|
if( |
| 553 |
|
$ (INDMAX(ic)-1).ge.INDSTART(ic) |
| 554 |
|
$ ) |
| 555 |
|
$ sl1 = CLSIGNAL(INDMAX(ic)-1) |
| 556 |
|
|
| 557 |
|
sl2 = -100000 !left 2 |
| 558 |
|
if( |
| 559 |
|
$ (INDMAX(ic)-2).ge.INDSTART(ic) |
| 560 |
|
$ ) |
| 561 |
|
$ sl2 = CLSIGNAL(INDMAX(ic)-2) |
| 562 |
|
|
| 563 |
|
sr1 = -100000 !right 1 |
| 564 |
|
if( |
| 565 |
|
$ (ic.ne.NCLSTR1.and.(INDMAX(ic)+1).lt.INDSTART(ic+1)) |
| 566 |
|
$ .or. |
| 567 |
|
$ (ic.eq.NCLSTR1.and.(INDMAX(ic)+1).le.TOTCLLENGTH) |
| 568 |
|
$ ) |
| 569 |
|
$ sr1 = CLSIGNAL(INDMAX(ic)+1) |
| 570 |
|
|
| 571 |
|
sr2 = -100000 !right 2 |
| 572 |
|
if( |
| 573 |
|
$ (ic.ne.NCLSTR1.and.(INDMAX(ic)+2).lt.INDSTART(ic+1)) |
| 574 |
|
$ .or. |
| 575 |
|
$ (ic.eq.NCLSTR1.and.(INDMAX(ic)+2).le.TOTCLLENGTH) |
| 576 |
|
$ ) |
| 577 |
|
$ sr2 = CLSIGNAL(INDMAX(ic)+2) |
| 578 |
|
|
| 579 |
|
if(ncog.ge.1)nbad = nbad+1-CLBAD(INDMAX(ic)) |
| 580 |
|
|
| 581 |
|
if(ncog.ge.2)then |
| 582 |
|
if(sl1.gt.sr1.and.(INDMAX(ic)-1).ge.1)then |
| 583 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)-1) |
| 584 |
|
elseif(sl1.le.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then |
| 585 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)+1) |
| 586 |
|
endif |
| 587 |
|
endif |
| 588 |
|
|
| 589 |
|
if(ncog.ge.3)then |
| 590 |
|
if(sl1.gt.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then |
| 591 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)+1) |
| 592 |
|
elseif(sl1.le.sr1.and.(INDMAX(ic)-1).ge.1)then |
| 593 |
|
c if(INDMAX(ic)-1.eq.0) |
| 594 |
|
c $ print*,' ======= ',sl2,sl1,sc,sr1,sr2 |
| 595 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)-1) |
| 596 |
|
endif |
| 597 |
|
endif |
| 598 |
|
|
| 599 |
|
if(ncog.ge.4)then |
| 600 |
|
if(sl2.gt.sr2.and.(INDMAX(ic)-2).ge.1)then |
| 601 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)-2) |
| 602 |
|
elseif(sl2.le.sr2.and.(INDMAX(ic)+2).le.NCLSTR1)then |
| 603 |
|
nbad=nbad+1-CLBAD(INDMAX(ic)+2) |
| 604 |
|
endif |
| 605 |
|
endif |
| 606 |
|
|
| 607 |
|
c if(ncog.ge.5)then |
| 608 |
|
c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG |
| 609 |
|
c $ ,' not implemented' |
| 610 |
|
c endif |
| 611 |
|
|
| 612 |
|
elseif(ncog.eq.0)then |
| 613 |
|
* ========================= |
| 614 |
|
* COG computation |
| 615 |
|
* ========================= |
| 616 |
|
|
| 617 |
|
|
| 618 |
|
|
| 619 |
|
iv=VIEW(ic) |
| 620 |
|
if(mod(iv,2).eq.1)incut=incuty |
| 621 |
|
if(mod(iv,2).eq.0)incut=incutx |
| 622 |
|
|
| 623 |
|
istart = INDSTART(IC) |
| 624 |
|
istop = TOTCLLENGTH |
| 625 |
|
if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1 |
| 626 |
|
nbad = 0 |
| 627 |
|
c$$$ do i=istart,istop |
| 628 |
|
c$$$ cut = incut*CLSIGMA(i) |
| 629 |
|
c$$$ if(CLSIGNAL(i).ge.cut)nbad = nbad +1 -CLBAD(i) |
| 630 |
|
c$$$ enddo |
| 631 |
|
do i = INDMAX(IC),istart,-1 |
| 632 |
|
cut = incut*CLSIGMA(i) |
| 633 |
|
if(CLSIGNAL(i).ge.cut)then |
| 634 |
|
nbad = nbad +1 -CLBAD(i) |
| 635 |
|
else |
| 636 |
|
goto 10 |
| 637 |
|
endif |
| 638 |
|
enddo |
| 639 |
|
10 continue |
| 640 |
|
do i = INDMAX(IC)+1,istop |
| 641 |
|
cut = incut*CLSIGMA(i) |
| 642 |
|
if(CLSIGNAL(i).ge.cut)then |
| 643 |
|
nbad = nbad +1 -CLBAD(i) |
| 644 |
|
else |
| 645 |
|
goto 20 |
| 646 |
|
endif |
| 647 |
|
enddo |
| 648 |
|
20 continue |
| 649 |
|
|
| 650 |
|
else |
| 651 |
|
|
| 652 |
|
c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG |
| 653 |
|
c $ ,' not implemented' |
| 654 |
|
|
| 655 |
|
|
| 656 |
|
endif |
| 657 |
|
|
| 658 |
|
nbadstrips = nbad |
| 659 |
|
|
| 660 |
|
return |
| 661 |
|
end |