chookEabsorb.f File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	chookeabsorbi (info)
subroutine	chookeabsorb (a, b, dE, info)
subroutine	chookeabsorbd (a, dE, info)
subroutine	chookeabsorbb (a, info)
subroutine	chookeabsorbc (a, n, p, info)
subroutine	chookeabsorbw (a, n, p, Eout, diff)

Function/Subroutine Documentation

chookeabsorb()

subroutine chookeabsorb	(	type(track)	a,
		type(track)	b,
		real*8	dE,
		integer	info
	)

Definition at line 51 of file chookEabsorb.f.

References wgt, and where.

      implicit none
#include "Zmaxdef.h"
#include "Zobs.h"
#include "Ztrack.h"
#include "Zabsorb.h"
!   This is called when Eabsorb != 0 and
!   when a charged particle runs from a 
!   to b and deposits energy dE (GeV) to the  Air.
!             
      type(track)::  a  ! input. charged particle track info. at a.
      type(track)::  b  ! input. charged particle track info. at b.
      real*8  de  ! input.  energy deposit in GeV; no weight is applied yet
      integer info ! in/out. not used now

      if(a.where .ge. 1)  then
         debydedx(a.where)= debydedx(a.where) + de*a.wgt
      endif

chookeabsorbb()

subroutine chookeabsorbb	(	type(track)	a,
		integer	info
	)

Definition at line 137 of file chookEabsorb.f.

References coszenith, false, mass, p, true, and wgt.

      implicit none
#include "Zmaxdef.h"
#include "Zcode.h"
#include "Ztrack.h"
#include "Ztrackp.h"
#include "Ztrackv.h"
#include "Zobs.h"
#include "Zobsp.h"
#include "Zobsv.h"
#include "Zabsorb.h"
!  This is called when Eabsorb(1) !=0 and
!    a particle crosses an observation level
!  info=2:  normal observation level
!      =1:  upper boundary (whcih is equal to or higher than the injection
!           height)
!      =3:  lower boundary (which is equal to or lower than the last
!           observation depth).  
!
      type(track)::  a  ! input.  a particle that croses the observation level
      integer info  ! input.  see above

      integer code
      integer lv
      logical count

      count =.false. 
      if(info .eq. 2 .and. a.where .eq. eabsorb(2) )  then
!           Eabsorb(2) is the layer number where the user want to
!           take sum of particle energy reaching there from above.
!                ptcl comes to the specified level OR
!                upper boundary.  This condition neglect
!                ptcles reaching the real lower boundary
!               (if Eabosrb(2) is not NoOfSites+1)
         if( a.vec.coszenith .gt. 0.)   then
!          if info=2,  this is to count energy  reaching
!          the level from above, we discard going up ptcls
!          (but it may come down, so some over count may happen) 
            count = .true.
            lv = 1
         endif
      elseif( info .eq. 1 ) then
         lv = 2      ! escape to space
         count = .true.
      endif
!
!            
      if(count) then
         code = a.p.code
!              for other ptcls than g,e,mu, pi, K, N, use 7
         if(code .gt. 7) code=7
!             at the last layer we see sum of the particle energy
!             what energy we should use here is somewhat annoying 
!             point.  We use total energy here.  
         if(lv .eq. 1) then
!            Ecrash(code) = Ecrash(code) +  a.p.fm.p(4)*a.wgt
            ecrash(code) = ecrash(code) + (a.p.fm.p(4)-a.p.mass)*a.wgt
         else
!           Espace(code) = Espace(code) + a.p.fm.p(4)*a.wgt
            espace(code) = espace(code) +  (a.p.fm.p(4)-a.p.mass)*a.wgt
         endif
      endif

chookeabsorbc()

subroutine chookeabsorbc	(	type(track)	a,
		integer	n,
		type(ptcl), dimension(n)	p,
		integer	info
	)

Definition at line 200 of file chookEabsorb.f.

References chookeabsorbw(), masn, and masp.

      use modxsecmedia
      implicit none
#include "Zmaxdef.h"
#include "Zcode.h"
#include "Zmass.h"
#include "Zobs.h"
#include "Ztrack.h"
#include "Ztrackv.h"
#include "Zabsorb.h"

      type(track)::  a    ! input.  incident particle track
                           !  that made the collision
      integer n     ! input. number of procuded ptcls in the collision
      type(ptcl)::  p(n) ! input. produced ptpcls.
      integer info  ! input. not used now.
      
      real*8 eout, diff, reldiff, ein
      real*8 diff1, diff2, ein1, ein2
      integer i
!          We check conservation above  this energy (GeV).
      real*8 ebig/5.d3/
      save
!          since target is not well known and  nucleon there
!          has Fermi momentum, we simply assume rest mass.
!          For E>5TeV, Ar mass ~40GeV/5000 < 4/500 < 1 % 
!          So we can see the conservation neglecting mass term
!            
      if( a.p.fm.p(4) .gt. ebig ) then
         ein1 = a.p.fm.p(4) + masn*(targetnucleonno-targetprotonno) +
     *       masp*targetprotonno
         ein2 = a.p.fm.p(4) + masp
         eout = 0.
         do i = 1, n
            eout = eout + p(i).fm.p(4)
         enddo

         diff1 = eout - ein1
         diff2 = eout - ein2
         ein = ein1
         diff = diff1
!              take smaller diff one for Ein
         if( abs(diff2) .lt. abs(diff1)) then
            diff = diff2
            ein = ein2
         endif
         reldiff = (eout/ein -1.)

         if( abs(reldiff) .gt. 0.02 ) then
            call chookeabsorbw(a, n, p, eout, diff)
         endif

         if( abs(diff) .gt. abs(maxebreak(1)) ) then
            maxebreak(1) = diff
            maxebreak(2) = reldiff
         endif

         if( abs(reldiff) .gt. abs(maxrelebreak(1)) ) then
            maxrelebreak(1) = reldiff
            maxrelebreak(2) = diff
         endif

         sumediff = sumediff + diff
         sumabsediff = sumabsediff + abs(diff)
      endif

Here is the call graph for this function:

chookeabsorbd()

subroutine chookeabsorbd	(	type(track)	a,
		real*8	dE,
		integer	info
	)

Definition at line 70 of file chookEabsorb.f.

References charge, code, kelec, kkaon, kneue, kneumu, knuc, kphoton, p, regptcl, subcode, wgt, and where.

      implicit none
#include "Zmaxdef.h"
#include "Zcode.h"
#include "Zobs.h"
#include "Ztrack.h"
#include "Zabsorb.h"
!   This is called when Eabsorb != 0 and
!   when a  particle energy becomes < Emin (info=0) for 
!   its traveling time becomes > limit or (info=2)
!   its angle relative to the 1ry becomes > limit.(info=4) 
!   Whether this is called or not depends on the   particle
!   and bit in Eabsorb.  dE is energy that can be regarded
!   as absorbed in the Air. (GeV) eventually.
!   bit 1 is the LSB of Eabsorb.
!     
!   bit   particle
!    1     photon: used to absorb shell energy at 
!          photoelectric effect. This bit is not used in the Air.  
!    2    photon.  
!    3    e+/e-
!    4    proton
!    5    neutron
!    6    anti-N
!    7    decaying prtcl 
!    8    others
!             
!***** Normally Eabsorb=6 (110 in bit pattarn) is enough.****
!
!
      type(track)::  a  ! input.  a particle that is < Emin 
                          ! at birth
      real*8  de  ! input.  energy which is supposed to  be emitted by 
                  ! the dying particle
      integer info ! in/out. not used now.
      if(a.where .ge. 1) then
         if(a.p.code .eq. kneue .or. a.p.code .eq. kneumu) then
           !  neutrino
            debydeathneu(a.where) = debydeathneu(a.where) + de*a.wgt
         elseif(a.p.code .eq. knuc .and. a.p.charge .eq. 0 .and.
     *          a.p.subcode .eq. regptcl ) then
            !  low E nutron
            debydeathnut(a.where) = debydeathnut(a.where) + de*a.wgt
         else
!                  next one and above are kept same as older versions for 
!                  compativilty 
            debydeath(a.where) = debydeath(a.where) + de*a.wgt

!                  we further put details 
            if(a.p.code .eq. kphoton ) then
               debydeathg(a.where) = debydeathg(a.where) + de*a.wgt
            elseif(a.p.code .eq. kelec ) then
               debydeathe(a.where) = debydeathe(a.where) + de*a.wgt
            elseif( a.p.code .le. kkaon )  then
               debydeathmupik(a.where) = debydeathmupik(a.where)
     *             + de*a.wgt
            elseif(a.p.code .eq. knuc .and. a.p.charge .eq. 1 ) then
!                 p
               debydeathp(a.where) = debydeathp(a.where) + de*a.wgt
            else
!                pbar, nbar, heavy,  others 
               debydeatho(a.where) = debydeatho(a.where) + de*a.wgt
            endif
         endif
      endif

chookeabsorbi()

subroutine chookeabsorbi

(

integer

info

)

Definition at line 2 of file chookEabsorb.f.

!           init for each event
      implicit none
#include "Zmaxdef.h"
#include "Zcode.h"
#include "Ztrack.h"
#include "Ztrackv.h"
#include "Zobs.h"
#include "Zobsp.h"
#include "Zobsv.h"
#include "Zabsorb.h"
      integer info ! not used now
      integer i
      do i = 0, noofsites + 1
!            i = 0 and NoOfSite+1 will not be used now. 
         debydedx(i) = 0.
         debydeath(i) = 0.
         debydeathg(i) = 0.
         debydeathe(i) = 0.
         debydeathmupik(i) = 0.
         debydeathneu(i) = 0.
         debydeathp(i) = 0.
         debydeathnut(i) = 0.
         debydeatho(i) = 0.
      enddo
      do i = 1, 7
!          i for g,e,mu,pi,K, n, (neu and others)
!             Energy crashing to the given layer
!             default is the NoOfSites. Can be changed by Eabosrb(2)
         ecrash(i) = 0.
!             Energy escaping to the upper bound( default is injection
!             height + 1m)  
         espace(i) = 0.
      enddo
!          counter used for checking energy conservation 
!          at the multiple producion.   
!           negative value is generated energy is < ininial energy
!           positive value is generated energy is > initial energy 
      do i = 1, 2
!         i=1 for max Ebreak, i=2 Relative Break value for that event
         maxebreak(i) = 0.
!         i=1 for max RelEgreak, i=2 E Break value for that event
         maxrelebreak(i) = 0.
      enddo

      sumediff = 0.
      sumabsediff = 0.

chookeabsorbw()

subroutine chookeabsorbw	(	type(track)	a,
		integer	n,
		type(ptcl), dimension(n)	p,
		real*8	Eout,
		real*8	diff
	)

Definition at line 266 of file chookEabsorb.f.

References code, cqeventno(), depth, height, subcode, and wgt.

      use modxsecmedia
      implicit none
#include "Zmaxdef.h"
#include "Zcode.h"
#include "Zmass.h"
#include "Zobs.h"
#include "Ztrack.h"
#include "Ztrackv.h"
#include "Zabsorb.h"

      type(track)::  a    ! input.  incident particle track
                           !  that made the collision
      integer n     ! input. number of procuded ptcls in the collision
      type(ptcl)::  p(n) ! input. produced ptpcls.
      real*8 eout, diff
      integer nevent, ntevent


      call cqeventno(nevent, ntevent)

      write(0,*) "=================  Warning: Ebreak ============="
      write(0,*) " event #=", nevent, " 1ry Energy=", inci.p.fm.p(4)
      write(0,*) "incident(not 1ry) code=", a.p.code,
     *     " subcode=", a.p.subcode, " Ein~", a.p.fm.p(4)
      write(0,*) "Target(A,Z)=", targetnucleonno, targetprotonno
      write(0,*) "No. of generted particles =", n
      write(0,*) "Sum of outgoing ptcl Energy:Eout=", eout
      write(0,*) "dE=(Eout-Ein)=", diff, " dE/Ein~",  
     *          diff/a.p.fm.p(4)
      write(0,*) "(Eout-Ein)/1ryE=",  diff/inci.p.fm.p(4)
      write(0,*) "Height=", a.pos.height," depth=",a.pos.depth/10.
      write(0,*) " where=", a.where, " weight=", a.wgt
      write(0,*) "================================================"

Here is the call graph for this function: