/[PAMELA software]/DarthVader/TrackerLevel2/src/F77/reductionflight.f
ViewVC logotype

Diff of /DarthVader/TrackerLevel2/src/F77/reductionflight.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.7 by pam-fi, Fri Sep 29 08:13:04 2006 UTC revision 1.19 by pam-fi, Mon May 14 11:03:06 2007 UTC
# Line 23  Line 23 
23        integer ierror        integer ierror
24        ierror = 0        ierror = 0
25    
26  *     -------------------------------------------------------  c$$$      debug = .true.
27  *     STRIP MASK  c$$$      verbose = .true.
28  *     -------------------------------------------------------  c$$$      warning = .true.
29    
30    *     //////////////////////////
31    *     initialize some parameters
32    *     //////////////////////////
33    
 c      call stripmask   !called later, after CN computation  
34        call init_level1        call init_level1
35    
36          if(debug)print*,'-- check LEVEL0 status'
37    
38  c      good1 = good0  c      good1 = good0
39  c--------------------------------------------------  c--------------------------------------------------
40  c     check the LEVEL0 event status for missing  c     check the LEVEL0 event status for missing
# Line 41  c--------------------------------------- Line 46  c---------------------------------------
46  c           ------------------------  c           ------------------------
47  c           GOOD  c           GOOD
48  c           ------------------------  c           ------------------------
49              GOOD1(DSPnumber(iv))=0                  !OK              GOOD1(DSPnumber(iv))=0 !OK
50  c           ------------------------  c           ------------------------
51  c           CRC error  c           CRC error
52  c           ------------------------  c           ------------------------
# Line 91  c           ------------------------ Line 96  c           ------------------------
96           good_old(iv)   = good1(iv)           good_old(iv)   = good1(iv)
97           ngood = ngood + good1(iv)           ngood = ngood + good1(iv)
98        enddo        enddo
99  c      if(ngood.ne.0)print*,'* WARNING * LEVEL0 event status: '        if(debug.and.ngood.ne.0)print*,'* WARNING * LEVEL0 event status: '
100  c     $     ,(good1(i),i=1,nviews)       $     ,(good1(i),i=1,nviews)
101  c--------------------------------------------------  c--------------------------------------------------
102  c     read the variable DATATRACKER from LEVEL0  c     read the variable DATATRACKER from LEVEL0
103  c     and fill the variable ADC (invertin view 11)  c     and fill the variable ADC (invertin view 11)
104  c--------------------------------------------------  c--------------------------------------------------
105          
106          if(debug)print*,'-- fill ADC vectors'
107    
108        call filladc(iflag)        call filladc(iflag)
109        if(iflag.ne.0)then        if(iflag.ne.0)then
 c        good1=0!<<<<<<<<<<<<<<<  
 c       if(DEBUG)print*,'event ',eventn(1),' >>>>>  decode ERROR'  
110           ierror = 220           ierror = 220
 c        goto 200  
 c         print*,'filladc error'  
111        endif        endif
112    
113  c--------------------------------------------------  c--------------------------------------------------
114  c     computes common noise for each VA1  c     computes common noise for each VA1
115  c     (excluding strips affected by signal,  c     (excluding strips with signal,
116  c     tagged with the flag CLSTR)  c     tagged with the flag CLSTR)
117  c--------------------------------------------------  c--------------------------------------------------
118          if(debug)print*,'-- compute CN'
119    
120        do iv=1,nviews        do iv=1,nviews
121          do ik=1,nva1_view           ima=0
122            cn(iv,ik)  = 0           do ik=1,nva1_view
123            cnn(iv,ik) = -1              cn(iv,ik)    = 0
124            mask_vk_ev(iv,ik)=1              cnrms(iv,ik) = 0
125            iflag=0              cnn(iv,ik)   = -1
126            if(mask_vk(iv,ik).eq.1)call cncomp(iv,ik,iflag)              iflag = 0
127  c          if(iflag.ne.0)good1=0              mask_vk_ev(iv,ik) = 1
128            if(iflag.ne.0)then              call stripmask(iv,ik)      !compute mask(i,j,k), combining VA1-masks
129               mask_vk_ev(iv,ik)=0  *           --------------------------------------
130               ierror = 220  *           if chip is not masked ---> evaluate CN
131               if(verbose)  *           --------------------------------------
132       $            print*,' * WARNING * Event ',eventn(1)              if( mask(iv,ik,1).eq.1 ) then !!!NBNB mask per la striscia 1 !!!!!!!!
133       $            ,': masked vk ',ik,' on view',iv                 call cncomp(iv,ik,iflag)
134            endif                 if(iflag.ne.0)then
135          enddo                    ima=ima+1
136                      mask_vk_ev(iv,ik)=0
137                      ierror = 220
138                   endif
139                   call stripmask(iv,ik) !compute mask(i,j,k), combining VA1-masks
140                endif
141             enddo
142     100     format(' * WARNING * Event ',i7,' view',i3,': VK MASK ',24i1)
143             if(ima.ne.0.and.debug)write(*,100)eventn(1),iv
144         $        ,(mask_vk_ev(iv,ik),ik=1,nva1_view)
145    c         if(ima.ne.0)write(*,100)eventn(1),iv
146    c     $        ,(mask_vk_ev(iv,ik),ik=1,nva1_view)  
147        enddo        enddo
 c      if(good1.eq.0)then  
 c         ierror = 220  
 c      endif  
148    
149        call stripmask !compute mask(i,j,k), combining mask_vk_ev and mask_vk  cc      call stripmask !compute mask(i,j,k), combining mask_vk_ev and mask_vk
150    
151  c---------------------------------------------  c---------------------------------------------
152  c     loops on views, VA1 and strips,  c     loops on views, VA1 and strips,
153  c     and computes strips signals using  c     and computes strips signals using
154  c     badstrip, pedestals, and  c     badstrip, pedestals, and
155  c     sigma informations from histograms  c     sigma informations from histograms
156  c---------------------------------------------  c---------------------------------------------
       flag_shower = .false.  
157        ind=1                     !clsignal array index        ind=1                     !clsignal array index
158    
159          if(debug)print*,'-- search clusters'
160        do iv=1,nviews            !loop on views        do iv=1,nviews            !loop on views
161          do is=1,nstrips_view    !loop on strips (1)          do is=1,nstrips_view    !loop on strips (1)
162            if(mod(iv,2).eq.1) then            if(mod(iv,2).eq.1) then
# Line 153  C===  > Y view Line 168  C===  > Y view
168       $           *mask(iv,nvk(is),nst(is))       $           *mask(iv,nvk(is),nst(is))
169              clinclcut(is)=incuty*sigma(iv,nvk(is),nst(is))              clinclcut(is)=incuty*sigma(iv,nvk(is),nst(is))
170       $           *mask(iv,nvk(is),nst(is))       $           *mask(iv,nvk(is),nst(is))
171  ccc            print*,"value(",is,")(reduction)= ",value(is)              sat(is)=0
172                if( adc(iv,nvk(is),nst(is)).lt.adc_saty )sat(is)=1
173            else                        else            
174  C===  > X view  C===  > X view
175              value(is)= (DBLE(adc(iv,nvk(is),nst(is)))              value(is)= (DBLE(adc(iv,nvk(is),nst(is)))
# Line 163  C===  > X view Line 179  C===  > X view
179       $           *mask(iv,nvk(is),nst(is))       $           *mask(iv,nvk(is),nst(is))
180              clinclcut(is)=incutx*sigma(iv,nvk(is),nst(is))              clinclcut(is)=incutx*sigma(iv,nvk(is),nst(is))
181       $           *mask(iv,nvk(is),nst(is))       $           *mask(iv,nvk(is),nst(is))
182                sat(is)=0
183                if( adc(iv,nvk(is),nst(is)).gt.adc_satx )sat(is)=1
184            endif            endif
 c$$$          print*,iv,is,' --- ',adc(iv,nvk(is),nst(is)),cn(iv,nvk(is))  
 c$$$     $         ,pedestal(iv,nvk(is),nst(is)),value(is)  
 c$$$     $         ,sigma(iv,nvk(is),nst(is))  
 c          if(value(is).gt.clseedcut(is))  
 c     $         print*,iv,is,' --- (ADC_PED_CN) ',value(is),clseedcut(is)  
185          enddo                   !end loop on strips (1)          enddo                   !end loop on strips (1)
186          call search_cluster(iv)          call search_cluster(iv)
187  c$$$        if(flag_shower.eqv..true.)then  
 c$$$          call init_level1                
 c$$$          good1=0  
 c$$$          goto 200              !jump to next event  
 c$$$        endif  
 ccc  
 ccc    modified by Elena (08/2006)  
 ccc  
188          if(.not.flag_shower)then          if(.not.flag_shower)then
189             call save_cluster(iv)             call save_cluster(iv)
190               if(debug)print*,'view ',iv,' #clusters ', nclstr_view
191          else          else
192             fshower(iv) = 1             fshower(iv) = 1
193             GOOD1(DSPn) = 11  c           GOOD1(DSPnumber(iv)) = 11 !AHAHAHAHA ORRORE!!
194               GOOD1(iv) = 11
195     101       format(' * WARNING * Event ',i7,' view',i3
196         $          ,' #clusters > ',i5,' --> MASKED')
197               if(debug)write(*,101)eventn(1),iv,nclstrmax_view
198          endif          endif
199        enddo                     ! end loop on views        enddo                     ! end loop on views
200        do iv=1,nviews        do iv=1,nviews
201          do ik=1,nva1_view          do ik=1,nva1_view
202            cnev(iv,ik)  = cn(iv,ik) !assigns computed CN to ntuple variables            cnev(iv,ik)    = cn(iv,ik) !assigns computed CN to ntuple variables
203            cnnev(iv,ik) = cnn(iv,ik) !assigns computed CN to ntuple variables            cnrmsev(iv,ik) = cnrms(iv,ik) !assigns computed CN to ntuple variables
204  ccc          print*,"cnev(",iv,",",ik,")(reduction)= ",cnev(iv,ik)            cnnev(iv,ik)   = cnn(iv,ik) !assigns computed CN to ntuple variables
205          enddo          enddo
206        enddo        enddo
207  C---------------------------------------------  C---------------------------------------------
# Line 203  C--------------------------------------- Line 214  C---------------------------------------
214        do iv = 1,nviews        do iv = 1,nviews
215           ngood = ngood + good1(iv)           ngood = ngood + good1(iv)
216        enddo        enddo
217        if(ngood.ne.0)print*,'* WARNING * Event ',eventn(1)        if(debug.and.ngood.ne.0)print*,'* WARNING * Event ',eventn(1)
218       $     ,':LEVEL1 event status: '       $     ,':LEVEL1 event status: '
219       $     ,(good1(i),i=1,nviews)       $     ,(good1(i),i=1,nviews)
220  c------------------------------------------------------------------------  c------------------------------------------------------------------------
# Line 246  c      good1 = 0 Line 257  c      good1 = 0
257           indmax(ic) = 0           indmax(ic) = 0
258           maxs(ic) = 0           maxs(ic) = 0
259           mult(ic) = 0                     mult(ic) = 0          
260           dedx(ic) = 0           sgnl(ic) = 0
261           whichtrack(ic) = 0           whichtrack(ic) = 0     !assigned @ level2
262    
263        enddo        enddo
264        do id=1,maxlength         !???        do id=1,maxlength         !???
# Line 267  c        crc1(iv)=0 Line 278  c        crc1(iv)=0
278                
279        return        return
280        end        end
281    
282  *---***---***---***---***---***---***---***---***  *---***---***---***---***---***---***---***---***
283  *  *
284  *  *
# Line 287  c        crc1(iv)=0 Line 299  c        crc1(iv)=0
299    
300  c     local variables  c     local variables
301        integer rmax,lmax         !estremi del cluster        integer rmax,lmax         !estremi del cluster
302        integer rstop,lstop       !per decidere quali strip includere nel cluster        integer rstop,lstop      
303                                  ! oltre il seed        integer first,last
304        integer first,last,diff   !per includere le strip giuste... !???        integer fsat,lsat
   
       integer multtemp          !temporary multiplicity variable  
305    
306        external nst        external nst
307    
 c------------------------------------------------------------------------  
 c     looks for clusters on each view  
 C     : CERCO STRIP SOPRA CLSEEDCUT, POI SCORRO A DX FINCHE'  
 c     NON TROVO  
 C     STRIP PIU' BASSA (in segnale/rumore)  
 C     => L'ULTIMA DELLA SERIE CRESCENTE  
 C     (LA PIU' ALTA) E' IL  
 C     CLUSTER SEED. POI SCORRO A SX E DX INCLUDENDO TUTTE  
 C     LE STRIP (FINO A 17 AL  
 C     MAX) CHE SUPERANO CLINCLCUT.  
 C     QUANDO CERCO IL CLUSTER SEED SUCCESSIVO SALTO LA STRIP  
 C     ADIACENTE A DESTRA  
 C     DELL'ULTIMO CLUSTER SEED (CHE SARA' NECESSARIAMENTE  
 C     PIU' BASSA) E PRENDO  
 C     COME SEED UNA STRIP SOLO SE IL SUO SEGNALE E'  
 C     MAGGIORE DI QUELLO DELLA STRIP  
 C     PRECEDENTE (PRATICAMENTE PER EVITARE CHE L'ULTIMA  
 C     STRIP DI UN GRUPPO DI STRIP  
 C     TUTTE SOPRA IL CLSEEDCUT VENGA AUTOMATICAMENTE PRESA  
 C     COME SEED... DEVE ESSERE  
 C     PRESA SOLO SE IL CLUSTER E' DOUBLE PEAKED...)  
 c------------------------------------------------------------------------  
 c     6 ottobre 2003  
 c     Elena: CLSEEDCUT = 7 (old value 10)  
 c     Elena: CLINCLCUT = 4 (old value 5)  
   
308        iseed=-999                !cluster seed index initialization        iseed=-999                !cluster seed index initialization
309    
310          inext=-999                !index where to start new cluster search
311    
312          flag_shower = .false.
313        nclstr_view=0        nclstr_view=0
314    
315        do jl=1,nladders_view     !1..3 !loops on ladders        do jl=1,nladders_view              !1..3 !loops on ladders
316           first=1+nstrips_ladder*(jl-1) !1,1025,2049  
317           last=nstrips_ladder*jl !1024,2048,3072           first = 1+nstrips_ladder*(jl-1) !1,1025,2049
318  c     X views have 1018 strips instead of 1024           last  = nstrips_ladder*jl       !1024,2048,3072
319    
320    *        X views have 1018 strips instead of 1024
321           if(mod(iv,2).eq.0) then           if(mod(iv,2).eq.0) then
322              first=first+3              first=first+3
323              last=last-3              last=last-3
# Line 336  c     X views have 1018 strips instead o Line 325  c     X views have 1018 strips instead o
325    
326           do is=first,last       !loop on strips in each ladder           do is=first,last       !loop on strips in each ladder
327    
328              if(is.le.iseed+1) goto 220  c---------------------------------------------
329  *******************************************************  c     new-cluster search starts at index inext
330  *     Elena 08/2006  c---------------------------------------------
331  * QUESTA PARTE NON E` ADEGUATA per cluster con grossi rilasci di carica              if(is.lt.inext) goto 220 ! next strip
 * perche` salva molte volte lo stesso cluster  
 * (salvo il cluster rispetto al primo massimo e basta...)  
 *******************************************************  
 c$$$c-----------------------------------------  
 c$$$c     after a cluster seed as been found,  
 c$$$c     look for next one skipping one strip on the right  
 c$$$c     (i.e. look for double peak cluster)  
 c$$$c-----------------------------------------  
 c$$$            if(is.ne.first) then  
 c$$$               if(value(is).le.value(is-1)) goto 220  
 c$$$            endif  
 c$$$c-----------------------------------------  
 c$$$c     skips cluster seed  
 c$$$c     finding if strips values are descreasing (a strip  
 c$$$c     can be a cluster seed only if previous strip value  
 c$$$c     is lower)  
 c$$$c-----------------------------------------  
 *******************************************************  
 * LA RICERCA PARTE DALL'ULTIMA STRIP SALVATA (***TEMPORANEO****)  
 *******************************************************  
             if(is.le.iseed+rmax+1) goto 220  
 *******************************************************  
332    
333              if(value(is).gt.clseedcut(is)) then              if(value(is).gt.clseedcut(is)) then
 ccc              print*,"value(",is,")=",value(is),  
 ccc     $             " .gt.clseedcut(",is,")=",clseedcut(is)  
334  c-----------------------------------------  c-----------------------------------------
335  c     possible SEED...  c     possible SEED...
336  c-----------------------------------------  c-----------------------------------------
337                 itemp=is                 itemp = is
338                   fsat = 0         ! first saturated strip
339                   lsat = 0         ! last saturated strip
340                 if(itemp.eq.last) goto 230 !estremo...                 if(itemp.eq.last) goto 230 !estremo...
341  ****************************************************  c              ------------------------                
342  *     modificato da Elena (08/2006) per salvare  c              search for first maximum
343  *     il cluster intorno al massimo assoluto  c              ------------------------
 ****************************************************  
 c$$$               do while(value(itemp)  
 c$$$     $              /sigma(iv,nvk(itemp),nst(itemp))  
 c$$$     $              .le.value(itemp+1)  
 c$$$     $              /sigma(iv,nvk(itemp+1),nst(itemp+1))) !BIAS: aggiustare il caso uguale!???  
 c$$$                  itemp=itemp+1  
 c$$$                  if(itemp.eq.last) goto 230 !stops if reaches last strip  
 c$$$               enddo            ! of the ladder  
344                 do while(                 do while(
345       $                   value(itemp).le.value(itemp+1)       $                   value(itemp).le.value(itemp+1)
346       $              .and.value(itemp+1).gt.clseedcut(itemp+1))       $              .and.value(itemp+1).gt.clseedcut(itemp+1))
347                    itemp=itemp+1                    itemp = itemp+1
348                    if(itemp.eq.last) goto 230 !stops if reaches last strip                    if(itemp.eq.last)   goto 230 !stops if reaches last strip
349                      if(sat(itemp).eq.1) goto 230 !stop if reaches a saturated strip
350                 enddo            ! of the ladder                 enddo            ! of the ladder
351   230           continue   230           continue
352  c-----------------------------------------  c              -----------------------------
353    c              check if strips are saturated
354    c              -----------------------------    
355                   if( sat(itemp).eq.1 )then
356                      fsat = itemp
357                      lsat = itemp
358                      if(itemp.eq.last) goto 231 !estremo...
359                      do while( sat(itemp+1).eq.1 )
360                         itemp = itemp+1
361                         lsat = itemp
362                         if(itemp.eq.last)   goto 231 !stops if reaches last strip
363                      enddo                  
364                   endif
365     231           continue
366    c---------------------------------------------------------------------------
367  c     fownd SEED!!!  c     fownd SEED!!!
368  c-----------------------------------------  c     (if there are saturated strips, the cluster is centered in the middle)
369                 iseed=itemp      c---------------------------------------------------------------------------
370  c----------------------------------------------------------                 if(fsat.eq.0.and.lsat.eq.0)then
371                      iseed = itemp ! <<< SEED
372                   else
373                      iseed = int((lsat+fsat)/2) ! <<< SEED
374    c$$$                  print*,'saturated strips ',fsat,lsat,iseed
375    c$$$                  print*,'--> ',(value(ii),ii=fsat,lsat)
376                   endif    
377    c---------------------------------------------------------------
378  c     after finding a cluster seed, checks also adjacent strips,  c     after finding a cluster seed, checks also adjacent strips,
379  C     and marks the ones exceeding clinclcut  C     and tags the ones exceeding clinclcut
380  c----------------------------------------------------------  c---------------------------------------------------------------
381                 ir=iseed         !indici destro                 ir=iseed         !indici destro
382                 il=iseed         ! e sinistro                 il=iseed         ! e sinistro
383                                
# Line 407  c--------------------------------------- Line 388  c---------------------------------------
388                 lstop=0          ! inclusion loop                 lstop=0          ! inclusion loop
389    
390                 do while(lstop.eq.0.or.rstop.eq.0) !shifts left and right from                 do while(lstop.eq.0.or.rstop.eq.0) !shifts left and right from
391                    ir=ir+1       !position index for strips on right side of  
392                                  ! cluster seed  
393                    il=il-1       !and for left side                    ir=ir+1       !index for right side
394                      il=il-1       !index for left side
395  c------------------------------------------------------------------------  c------------------------------------------------------------------------
396  c     checks for last or first strip of the ladder  c     checks for last or first strip of the ladder
397  c------------------------------------------------------------------------  c------------------------------------------------------------------------
398                    if(ir.gt.last) then !when index goes beyond last strip                    if( ir.gt.last  ) rstop = 1                      
399                       rstop=1    ! of the ladder, change rstop flag in order                    if( il.lt.first ) lstop = 1
                                 ! to "help" exiting from loop  
                   endif  
400                                        
                   if(il.lt.first) then !idem when index goes beyond  
                      lstop=1    ! first strip of the ladder  
                   endif  
                     
 c------------------------------------------------------------------------  
 c     check for clusters including more than nclstrp strips  
401  c------------------------------------------------------------------------  c------------------------------------------------------------------------
402                    if((rmax-lmax+1).ge.nclstrp) then  c     add strips exceeding inclusion cut
                      goto 210   !exits inclusion loop:  
                                 ! lmax and rmax maintain last value  
                                 ! NB .ge.!???  
                   endif  
 c------------------------------------------------------------------------  
 c     marks strips exceeding inclusion cut  
403  c------------------------------------------------------------------------  c------------------------------------------------------------------------
404                    if(rstop.eq.0) then !if last strip of the ladder or last                    if( (rmax-lmax+1).ge.nclstrp )goto 210   !exits inclusion loop
405                                  ! over-cut strip has not been reached  
406                       if(value(ir).gt.clinclcut(ir)) then !puts in rmax the                    if(rstop.eq.0) then !if right cluster morder has not been reached
407                          rmax=ir ! last right over-cut strip                       if(value(ir).gt.clinclcut(ir)) then
408                            rmax=ir !include a strip on the right
409                       else                       else
410                          rstop=1 !otherwise cluster ends on right and rstop                          rstop=1 !cluster right end
411                       endif      ! flag=1 signals it                       endif    
412                    endif                    endif
413                    if(lstop.eq.0) then  
414                      if( (rmax-lmax+1).ge.nclstrp )goto 210   !exits inclusion loop
415    
416                      if(lstop.eq.0) then !if left cluster morder has not been reached
417                       if(value(il).gt.clinclcut(il)) then                       if(value(il).gt.clinclcut(il)) then
418                          lmax=il                          lmax=il !include a strip on the left
419                       else                       else
420                          lstop=1                          lstop=1 !cluster left end
421                       endif                       endif
422                    endif                    endif
423    
424                 enddo            !ends strip inclusion loop                 enddo            !ends strip inclusion loop
425                   goto 211
426   210           continue         !jumps here if more than nclstrp have been included   210           continue         !jumps here if more than nclstrp have been included
427                        c               print*,'>>> nclstrp! '
428                 multtemp=rmax-lmax+1 !stores multiplicity in temp   211           continue
429                                  ! variable. NB rmax and lmax can change later in  c-----------------------------------------
430                                  ! order to include enough strips to calculate eta3  c     end of inclusion loop!
431                                  ! and eta4. so mult is not always equal to cllength  c-----------------------------------------
432  c------------------------------------------------------------------------                
433  c     NB per essere sicuro di poter calcolare eta3 e eta4 devo includere  c------------------------------------------------------------------------
434  c     sempre e comunque le 2 strip adiacenti al cluster seed e quella  c     adjust the cluster in order to have at least a strip around the seed
435  c     adiacente ulteriore dalla parte della piu' alta fra queste due  c------------------------------------------------------------------------
436  c     (vedi oltre...)!???                 if(iseed.eq.lmax.and.lmax.ne.first)then
437  c------------------------------------------------------------------------                    lmax = lmax-1
438                      if( (rmax-lmax+1).gt.nclstrp )rmax=rmax-1
439  c     nel caso di estremi del ladder...!???                 endif
440                   if(iseed.eq.rmax.and.rmax.ne.last )then
441  c     ho meno di 4 strip nel cluster --> se sono sui bordi o quasi del ladder                    rmax = rmax+1
442  c     costruisco il cluster ad hoc e poi esco, se non sono sui bordi o quasi                    if( (rmax-lmax+1).gt.nclstrp )lmax=lmax+1
 c     vado oltre (aggiungero' quindi strip a sx e dx in modo da poter calcolare  
 c     eta3e4)  
                if((rmax-lmax+1).lt.4) then  
   
                   if(iseed.eq.first) then !estremi...  
                      rmax=iseed+2 !NB in questo modo puo' anche capitare di  
                      lmax=iseed ! includere strip sotto taglio di inclusione  
                      goto 250   ! che non serviranno per eta3e4!???  
                   endif  
                     
                   if(iseed.eq.last) then !estremi...  
                      rmax=iseed  
                      lmax=iseed-2 !NB 2 e non 3, perche' altrimenti sarei in  
                      goto 250   ! ((rmax-lmax+1).lt.4).eq.false. !???  
                   endif         !NMB questo e' l'unico caso di cllength=3!???  
                     
                   if(iseed.eq.first+1) then !quasi estremi...  
                      rmax=iseed+2  
                      lmax=iseed-1  
                      goto 250  
                   endif  
                   if(iseed.eq.last-1) then  
                      rmax=iseed+1  
                      lmax=iseed-2  
                      goto 250  
                   endif  
 c     se ho 4 o piu' strip --> se sono sui bordi esco, se sono sui quasi bordi  
 c     includo la strip del bordo  
                else  
                     
                   if(iseed.eq.first) goto 250 !estremi... non includo altro                    
                   if(iseed.eq.last) goto 250  
                   if(iseed.eq.first+1) then !quasi estremi... mi assicuro di  
                      lmax=first ! avere le strip adiacenti al seed  
                      if((rmax-lmax+1).gt.nclstrp) rmax=rmax-1 !NB effetto  
                      goto 250   ! coperta: se la lunghezza del cluster era gia'  
                   endif         ! al limite (nclstrp), per poter aggiungere questa  
                                 ! strip a sinistra devo toglierne una a destra...!???  
                   if(iseed.eq.last-1) then  
                      rmax=last  
                      if((rmax-lmax+1).gt.nclstrp) lmax=lmax+1  
                      goto 250  
                   endif                    
443                 endif                 endif
444    
445  c------------------------------------------------------------------------  c------------------------------------------------------------------------
446  c     be sure to include in the cluster the cluster seed with its 2 adjacent  c     adjust the cluster in order to store a minimum number of strips
447  c     strips, and the one adjacent to the greatest between this two strip, as the  c------------------------------------------------------------------------
448  c     fourth one. if the strips have the same value (!) the fourth one is chosen                 do while( (rmax-lmax+1).lt.nclstrpmin )
449  c     as the one having the greatest value between the second neighbors  
450  c------------------------------------------------------------------------                    vl = -99999
451                 if(value(iseed+1).eq.value(iseed-1)) then                    vr = -99999
452                    if(value(iseed+2).ge.value(iseed-2)) then !??? qui cmq c'e'                    if(lmax-1.ge.first) vl = value(lmax-1)
453                       diff=(iseed+2)-rmax                    if(rmax+1.le.last ) vr = value(rmax+1)
454                       if(diff.gt.0) then                    if(vr.ge.vl) then
455                          rmax=rmax+diff                       rmax = rmax+1
456                          if((rmax-lmax+1).gt.nclstrp) then                    else  
457                             lmax=rmax-nclstrp+1                       lmax = lmax-1
                         endif  
                      endif  
                      diff=(iseed-1)-lmax  
                      if(diff.lt.0) then  
                         lmax=lmax+diff  
                         if((rmax-lmax+1).gt.nclstrp) then  
                            rmax=lmax+nclstrp-1  
                         endif  
                      endif  
                   else  
                      diff=(iseed-2)-lmax  
                      if(diff.lt.0) then  
                         lmax=lmax+diff  
                         if((rmax-lmax+1).gt.nclstrp) then  
                            rmax=lmax+nclstrp-1  
                         endif  
                      endif  
                      diff=(iseed+1)-rmax  
                      if(diff.gt.0) then  
                         rmax=rmax+diff  
                         if((rmax-lmax+1).gt.nclstrp) then  
                            lmax=rmax-nclstrp+1  
                         endif  
                      endif  
                   endif  
                elseif(value(iseed+1).gt.value(iseed-1)) then  
 c     !??? sposto il limite del cluster a destra per includere sempre le strip  
 c     necessarie al calcolo di eta-i  
 c     se il cluster diventa  troppo lungo lo accorcio a sinistra per avere non piu'  
 c     di nclstrp (in questo caso sono sicuro di aver gia' incluso le strip  
 c     necessarie al calcolo di eta-i a sinistra, quindi se voglio posso uscire)  
                   diff=(iseed+2)-rmax  
                   if(diff.gt.0) then  
                      rmax=rmax+diff  
                      if((rmax-lmax+1).gt.nclstrp) then  
                         lmax=rmax-nclstrp+1  
 c     goto 250  
                      endif  
                   endif  
                   diff=(iseed-1)-lmax  
                   if(diff.lt.0) then  
                      lmax=lmax+diff  
                      if((rmax-lmax+1).gt.nclstrp) then  
                         rmax=lmax+nclstrp-1  
 c     goto 250 !inutile!???  
                      endif  
                   endif  
                else  
                   diff=(iseed-2)-lmax  
                   if(diff.lt.0) then  
                      lmax=lmax+diff  
                      if((rmax-lmax+1).gt.nclstrp) then  
                         rmax=lmax+nclstrp-1  
 c     goto 250  
                      endif  
                   endif  
                   diff=(iseed+1)-rmax  
                   if(diff.gt.0) then  
                      rmax=rmax+diff  
                      if((rmax-lmax+1).gt.nclstrp) then  
                         lmax=rmax-nclstrp+1  
 c     goto 250 !inutile!???  
                      endif  
458                    endif                    endif
459                 endif                    
460   250           continue                 enddo
461    
462  c--------------------------------------------------------  c--------------------------------------------------------
463  c     fills cluster variables  c     store cluster info
464  c--------------------------------------------------------  c--------------------------------------------------------
 c$$$               nclstr1=nclstr1+1 !cluster number  
 c$$$ccc               print*,nclstr1,multtemp  
 c$$$               if(nclstr1.gt.nclstrmax) then !too many clusters for the event:  
 c$$$                  if(verbose)print*,'Event ',eventn(1),  
 c$$$     $                 ': more than ',nclstrmax,' clusters'  
 c$$$                  good1=0       ! event  
 c$$$                  nclstr1=0  
 c$$$                  totCLlength=0  
 c$$$                  flag_shower = .true.  
 c$$$                  goto 2000  
 c$$$               endif  
 c$$$               view(nclstr1)   = iv !vista del cluster  
 c$$$               ladder(nclstr1) = nld(iseed,iv) !ladder a cui appartiene il cluster seed  
 c$$$               maxs(nclstr1)   = iseed !strip del cluster seed  
 c$$$               mult(nclstr1)   = multtemp !molteplicita'  
 c$$$                
 c$$$               indstart(nclstr1) = ind !posizione dell'inizio del cluster nell'  
 c$$$c                                      ! array clsignal  
 c$$$               indmax(nclstr1)   = indstart(nclstr1)+(iseed-lmax) !posizione del  
 c$$$c                                      ! cluster seed nell'array clsignal  
 c$$$                
 c$$$               CLlength      = rmax-lmax+1 !numero di strip del cluster  
 c$$$               totCLlength   = totCLlength+CLlength  
 c$$$               dedx(nclstr1) = 0  
 c$$$               do j=lmax,rmax   !stores sequentially cluter strip values in  
 c$$$                  clsignal(ind) = value(j) ! clsignal array  
 c$$$                  ind=ind+1  
 c$$$c                  if(value(j).gt.0)  
 c$$$                  if(value(j).gt.clinclcut(j))  
 c$$$     $                 dedx(nclstr1) = dedx(nclstr1) + value(j) !cluster charge  
 c$$$               enddo  
 ccc  
 ccc            *** Modified by Elena (08/2006) ***  
 ccc  
465                 nclstr_view = nclstr_view + 1 !cluster number                 nclstr_view = nclstr_view + 1 !cluster number
466  c               print*,'view ',iv,' -- search_cluster -- nclstr_view: '  
 c     $              ,nclstr_view  
467                 if(nclstr_view.gt.nclstrmax_view) then !too many clusters for the view:                 if(nclstr_view.gt.nclstrmax_view) then !too many clusters for the view:
468                    if(verbose) print*,'Event ',eventn(1),  c$$$                  if(verbose) print*,'Event ',eventn(1),
469       $                 ': more than ',nclstrmax_view  c$$$     $                 ': more than ',nclstrmax_view
470       $                 ,' clusters on view ',iv  c$$$     $                 ,' clusters on view ',iv
 c                  good1=0       ! event  
 c                  nclstr1=0  
 c                  totCLlength=0  
471                    flag_shower = .true.                    flag_shower = .true.
472                    goto 2000                    goto 2000
473                 endif                 endif
474    
475  c               view(nclstr1)   = iv !vista del cluster                 ladder_view(nclstr_view) = nld(iseed,iv)
476                 ladder_view(nclstr_view) = nld(iseed,iv) !ladder a cui appartiene il cluster seed                 maxs_view(nclstr_view)   = iseed
477                 maxs_view(nclstr_view)   = iseed !strip del cluster seed                 mult_view(nclstr_view)   = rmax-lmax+1
                mult_view(nclstr_view)   = multtemp !molteplicita'  
478                 rmax_view(nclstr_view)   = rmax                 rmax_view(nclstr_view)   = rmax
479                 lmax_view(nclstr_view)   = lmax                 lmax_view(nclstr_view)   = lmax
480    
481    c$$$               if(rmax-lmax+1.gt.25)
482    c$$$     $              print*,'view ',iv
483    c$$$     $              ,' cl ',nclstr_view,' mult ',rmax-lmax+1
484    c------------------------------------------------------------------------
485    c     search for a double peak inside the cluster                                                                                                            
486    c------------------------------------------------------------------------
487                   inext = rmax+1   !<< index where to start new-cluster search
488                  
489                   vmax = 0
490                   vmin = value(iseed)
491                   imax = iseed
492                   imin = iseed
493                   do iss = max(iseed+1,lsat+1),rmax
494                      if( value(iss).lt.vmin )then
495                         if( imax.ne.iseed )goto 221 !found dowble peek
496                         imin = iss
497                         vmin = value(iss)
498                      else
499                         delta = value(iss) - value(imin)
500                         cut = sqrt(clinclcut(iss)**2 + clinclcut(imin)**2)
501                         if(
502         $                    delta.gt.cut .and.
503         $                    value(iss).gt.clseedcut(iss).and.
504         $                    .true.)then
505                            if( value(iss).gt.vmax )then                        
506                               imax = iss
507                               vmax = value(iss)
508                            else
509                               goto 221 !found dowble peek
510                            endif
511                         endif
512                      endif
513                   enddo
514     221           continue
515                  
516                   if(imax.gt.iseed)then
517                      inext = imax    !<< index where to start new-cluster search
518    c$$$                  print*,'--- double peek ---'
519    c$$$                  print*,(value(ii),ii=lmax,rmax)
520    c$$$                  print*,'seed ',iseed,' imin ',imin,' imax ',imax
521                   endif
522  c--------------------------------------------------------  c--------------------------------------------------------
523  c  c
524  c--------------------------------------------------------  c--------------------------------------------------------
# Line 698  c        posizione del cluster seed nell Line 568  c        posizione del cluster seed nell
568                    
569           CLlength      = rmax_view(ic) - lmax_view(ic) + 1 !numero di strip salvate           CLlength      = rmax_view(ic) - lmax_view(ic) + 1 !numero di strip salvate
570           totCLlength   = totCLlength + CLlength           totCLlength   = totCLlength + CLlength
571           dedx(nclstr1) = 0           sgnl(nclstr1) = 0
572           do j=lmax_view(ic),rmax_view(ic)         !stores sequentially cluter strip values in           do j=lmax_view(ic),rmax_view(ic)         !stores sequentially cluter strip values in
573    
574              clsignal(ind) = value(j) ! clsignal array              clsignal(ind) = value(j) ! clsignal array
# Line 714  c            clped(ind)   = pedestal(iv, Line 584  c            clped(ind)   = pedestal(iv,
584              ind=ind+1              ind=ind+1
585  c     if(value(j).gt.0)  c     if(value(j).gt.0)
586              if(value(j).gt.clinclcut(j))              if(value(j).gt.clinclcut(j))
587       $           dedx(nclstr1) = dedx(nclstr1) + value(j) !cluster charge       $           sgnl(nclstr1) = sgnl(nclstr1) + value(j) !cluster charge
588           enddo           enddo
589    
590  c         print*,'view ',iv,' -- save_cluster -- nclstr1: '  c         print*,'view ',iv,' -- save_cluster -- nclstr1: '
591  c     $        ,nclstr1,maxs(nclstr1),mult(nclstr1),dedx(nclstr1)  c     $        ,nclstr1,maxs(nclstr1),mult(nclstr1),sgnl(nclstr1)
592                    
593        enddo        enddo
594                
# Line 733  c     $        ,nclstr1,maxs(nclstr1),mu Line 603  c     $        ,nclstr1,maxs(nclstr1),mu
603  *---***---***---***---***---***---***---***---***  *---***---***---***---***---***---***---***---***
604    
605    
606        subroutine stripmask  c$$$      subroutine stripmask
607    c$$$
608    c$$$*     this routine set va1 and single-strip masks,
609    c$$$*     on the basis of the VA1 mask saved in the DB
610    c$$$*
611    c$$$*     mask(nviews,nva1_view,nstrips_va1) !strip mask
612    c$$$*     mask_vk(nviews,nva1_view)          !VA1 mask
613    c$$$*
614    c$$$      include 'commontracker.f'
615    c$$$      include 'level1.f'
616    c$$$      include 'common_reduction.f'
617    c$$$      include 'calib.f'
618    c$$$
619    c$$$*     init mask
620    c$$$      do iv=1,nviews
621    c$$$         do ivk=1,nva1_view
622    c$$$            do is=1,nstrips_va1
623    c$$$c               mask(iv,ivk,is) = mask_vk(iv,ivk)
624    c$$$               if( mask_vk(iv,ivk) .ne. -1)then
625    c$$$                  mask(iv,ivk,is) = 1
626    c$$$     $                 * mask_vk(iv,ivk)     !from DB
627    c$$$     $                 * mask_vk_ev(iv,ivk)  !from <SIG>
628    c$$$     $                 * mask_vk_run(iv,ivk) !from CN
629    c$$$               else
630    c$$$                  mask(iv,ivk,is) = -1
631    c$$$     $                 * mask_vk(iv,ivk)     !from DB
632    c$$$     $                 * mask_vk_ev(iv,ivk)  !from CN
633    c$$$               endif
634    c$$$            enddo
635    c$$$         enddo
636    c$$$      enddo
637    c$$$
638    c$$$
639    c$$$      return
640    c$$$      end
641    
642          subroutine stripmask(iv,ivk)
643    
644    *     -----------------------------------------------
645  *     this routine set va1 and single-strip masks,  *     this routine set va1 and single-strip masks,
646  *     on the basis of the VA1 mask saved in the DB  *     on the basis of the VA1 mask saved in the DB
647  *  *
648  *     mask(nviews,nva1_view,nstrips_va1) !strip mask  *     mask(nviews,nva1_view,nstrips_va1) !strip mask
649  *     mask_vk(nviews,nva1_view)          !VA1 mask  *     mask_vk(nviews,nva1_view)          !VA1 mask
650  *  *     -----------------------------------------------
651        include 'commontracker.f'        include 'commontracker.f'
652        include 'level1.f'        include 'level1.f'
653        include 'common_reduction.f'        include 'common_reduction.f'
654        include 'calib.f'        include 'calib.f'
655    
656  *     init mask  *     init mask
657        do iv=1,nviews        do is=1,nstrips_va1
658           do ivk=1,nva1_view  *        --------------------------------------------------------
659              do is=1,nstrips_va1  *        if VA1-mask from DB is 0 or 1, three masks are combined:
660  c               mask(iv,ivk,is) = mask_vk(iv,ivk)  *        - from DB (a-priori mask)
661                 mask(iv,ivk,is) = mask_vk(iv,ivk) * mask_vk_ev(iv,ivk)                *        - run-based (chip declared bad on the basis of <SIG>)
662              enddo  *        - event-based (failure in CN computation)
663           enddo  *        --------------------------------------------------------
664             if( mask_vk(iv,ivk) .ne. -1)then            
665                mask(iv,ivk,is) = 1
666         $           * mask_vk(iv,ivk)     !from DB
667         $           * mask_vk_ev(iv,ivk)  !from <SIG>
668         $           * mask_vk_run(iv,ivk) !from CN
669    *        -----------------------------------------------------------
670    *        if VA1-mask from DB is -1 only event-based mask is applied
671    *        -----------------------------------------------------------
672             else
673                mask(iv,ivk,is) = -1
674         $           * mask_vk(iv,ivk)     !from DB
675         $           * mask_vk_ev(iv,ivk)  !from CN
676             endif
677        enddo        enddo
678          
679          
680        return        return
681        end        end
   

Legend:
Removed from v.1.7  
changed lines
  Added in v.1.19

  ViewVC Help
Powered by ViewVC 1.1.23