crigCut.f File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	crigcut (azmin, zen, rig, prob)
subroutine	k2dtblv (tbl, xs, ys, adj, xm, ym, dx, dy, x, y, ans)
subroutine	cgetrigcut (azmin, zen, rigcut)
subroutine	cgetrigprob (tbl, izen, iphi, irig, i1, i2, i3, prob)
subroutine	crigcut0 (file)
subroutine	crigcutfmt1
subroutine	creadrigcut (io, tbl, azm, zen, adj)
subroutine	cfillrigcut (tbl, azm, zen)
subroutine	cqrigcutplace (tlt, tlg, mdec)
subroutine	crigcutfmt2
subroutine	crigcutfmt4
subroutine	crdrigcut2 (fmt, io, tbl, izen, iphi, irig)
subroutine	crdrigcut4 (fmt, io, minval, tbl, izen, irig)
subroutine	crigcutfmt5

Function/Subroutine Documentation

cfillrigcut()

subroutine cfillrigcut	(	real*8, dimension(azm, zen)	tbl,
		integer	azm,
		integer	zen
	)

Definition at line 284 of file crigCut.f.

Referenced by crigcutfmt1(), and crigcutfmt5().

       implicit none
       integer  azm, zen
       real*8 tbl(azm,  zen)
       integer i

        do i = 1, zen
            tbl(azm, i) = tbl(1, i)
        enddo

Here is the caller graph for this function:

cgetrigcut()

subroutine cgetrigcut	(	real*8	azmin,
		real*8	zen,
		real*8	rigcut
	)

Definition at line 117 of file crigCut.f.

References cerrormsg(), d0, and k4ptdi().

Referenced by __drawrigcut.f__().

      implicit none
#include "ZrigCut.h"
!         get rigidity cut (for old table only); use crigCut now
!
      real*8 azmin  ! input. value of azimuth angle at a given location.
                   !         deg.
      real*8 zen   ! input. value of zenith angle at a given location. 
      real*8 rigcut ! output. cutoff rigidity. in GV
!

      real*8 azm

      azm = azmin

      if(azm .lt. 0.) then
         azm = azm + 360.
      endif
      azm = mod(azmin, 360.d0)
      if(rdatafmt .eq. 1) then
         call k4ptdi(rigcuttbl, azmsize, zensize, azmsize, minazm, 
     *        minzen, dazm,  dzen,  azm, zen, rigcut)
      else
         call cerrormsg(
     *       'only old cutoff table can be used', 0)
      endif

Here is the call graph for this function:

Here is the caller graph for this function:

cgetrigprob()

subroutine cgetrigprob	(	real, dimension(izen, iphi, irig)	tbl,
		integer	izen,
		integer	iphi,
		integer	irig,
		integer	i1,
		integer	i2,
		integer	i3,
		real*8	prob
	)

Definition at line 145 of file crigCut.f.

Referenced by crigcut().

      implicit none
#include "ZrigCut.h"

      integer izen, iphi, irig
      real tbl(izen, iphi, irig)
      integer i1, i2, i3
      real*8 prob

!      if( i3 .lt. RigSize ) then
!       prob = dlogrig/LogDRig * 
!     *          (tbl(i1,i2,i3+1)-tbl(i1, i2, i3))
!     *         +  tbl(i1,i2,i3)
!      else
!         prob = tbl(i1,i2,i3)
!      endif

      prob = tbl(i1,i2,i3)

Here is the caller graph for this function:

cqrigcutplace()

subroutine cqrigcutplace	(	real*8	tlt,
		real*8	tlg,
		real*8	mdec
	)

Definition at line 295 of file crigCut.f.

       implicit none
#include "ZrigCut.h"
       real*8 tlt, tlg, mdec


!            to inform lat, long, magdec
        tlt = latit
        tlg = longi
        mdec = magdec

crdrigcut2()

subroutine crdrigcut2	(	integer	fmt,
		integer	io,
		real*4, dimension(izen, iphi, irig)	tbl,
		integer	izen,
		integer	iphi,
		integer	irig
	)

Definition at line 396 of file crigCut.f.

References cerrormsg().

Referenced by crigcutfmt2().

!       read cut off table for format 2
      implicit none
      integer fmt ! input.  fmt=2 or 3, if  3, data for izen, iphi, irig are
                  !          not given
      integer io  ! input. table dev. number
      integer izen ! input. table for zenith angle
      integer iphi ! input. table for azimuthal angle
      integer irig ! input. table for rigidity
      real*4  tbl(izen, iphi, irig)  ! output. 3D table which shows prob. that
                          ! given ptcl with izen, iphi, irig can enter the earth
                  !  NOTE:  this  is  real*4--------------------

      integer i1, i2, i3
      integer j1, j2, j3
      character*128 msg


      do i1= 1 , izen
         do i2 = 1, iphi
            do i3 = 1, irig
               if(fmt .eq. 2) then
                  read(io, *, end=100) j1, j2, j3, tbl(i1, i2, i3)
                  if(i1 .ne. j1+1 .or.
     *               i2 .ne. j2+1 .or.
     *                       i3 .ne. j3+1) then 
                     write(msg, *)
     *               ' data index mismatch in new rigidit cut table',
     *                 i1, j1+1, i2, j2+1, i2, j3+1      
                     call cerrormsg(msg, 0)
                  endif
               else
                  read(io, *) tbl(i1, i2, i3)
               endif
            enddo
         enddo
      enddo
      return
 100  continue
      call cerrormsg(
     *   'Unexpected EOF in new rigidity cut table',0)

Here is the call graph for this function:

Here is the caller graph for this function:

crdrigcut4()

subroutine crdrigcut4	(	integer	fmt,
		integer	io,
		real*8	minval,
		real*4, dimension(izen, irig)	tbl,
		integer	izen,
		integer	irig
	)

Definition at line 439 of file crigCut.f.

References cerrormsg(), and e.

Referenced by crigcutfmt4().

!       read cut off table for format 4
      implicit none
      integer fmt ! input.  fmt=4
      integer io  ! input. table dev. number
      real*8  minval ! input minimum value of rigidty (for check)
      integer izen ! input. table for zenith angle
      integer irig ! input. table for rigidity
      real*4  tbl(izen, irig)  ! output. 2D table which shows prob. that
                          ! given ptcl with izen, irig can enter the earth
                  !  NOTE:  this  is  real*4--------------------

      integer i1, i3
      integer j1,  idummy
      real*4 dummy
      character*128 msg


      do i1= 1 , izen
         do i3 = 1, irig
            read(io, *, end=100)
     *         j1, idummy, idummy, dummy, tbl(i1, i3)
            if(i3 .eq. 1) then
               if( abs(dummy - minval)/minval .gt. 1.e-3) then
                  write(msg,*)
     *             'check min rigidity in headr=', minval,
     *             ' in table=', dummy
                  call cerrormsg(msg, 0)
               endif
            endif
            if(i1 .ne. j1) then
               write(msg, *)
     *           ' data index mismatch in new rigidit cut table',
     *              i1, j1
               call cerrormsg(msg, 0)
            endif
         enddo
      enddo
      return
 100  continue
      call cerrormsg(
     *   'Unexpected EOF in new rigidity cut table',0)

Here is the call graph for this function:

Here is the caller graph for this function:

creadrigcut()

subroutine creadrigcut	(	integer	io,
		real*8, dimension(adj, zen)	tbl,
		integer	azm,
		integer	zen,
		integer	adj
	)

Definition at line 263 of file crigCut.f.

References cerrormsg().

Referenced by crigcutfmt1(), and crigcutfmt5().

       implicit none
       integer io,  azm, zen, adj
       real*8 tbl(adj,  zen)


       integer i, j, ios
       character*150 msg

       do   j=1,  zen
           read(io, *, iostat=ios) ( tbl(i, j), i = 1, azm)
           if(ios .ne. 0) then
              write(msg, *) 'Unexpected EOF at rigicity table reading'
              call cerrormsg(msg, 1)
              write(msg, *) ' line number=', j, ' azm=',azm, ' zen=',
     *          zen, ' adj=',adj
              call cerrormsg(msg, 0)
           endif
       enddo

Here is the call graph for this function:

Here is the caller graph for this function:

crigcut()

subroutine crigcut	(	real*8	azmin,
		real*8	zen,
		real*8	rig,
		real*8	prob
	)

Definition at line 6 of file crigCut.f.

References cerrormsg(), cgetrigprob(), d0, k2dtblv(), and k4ptdi().

Referenced by cifcutoff(), crigfunc2(), funcazim(), and primdn().

      implicit none
#include "ZrigCut.h"
!          Before calling this subroutine, you must have called
!       crigCut0 to get the table which contains the information
!       of cutoff rigidities.

!
      real*8 azmin  ! input. value of azimuth angle at a given location.
                   !         deg.
                   ! For format# 5, this is longitude of a place
      real*8 zen   ! input. value of zenith angle at a given location. 
                   !     It may be deg, or cos, depending on the table.
                   ! For format# 5, this is latitude of the place.
                   !     It may be deg, or cos, depending on the table.
      real*8 rig   ! input. value of particle rigidity (GV) p/Z

      real*8 prob  ! output. probability that 'r' can enter the earth.
                   !  For format# 5, this is  0 or 1. 1 dose not necessary
                   !   mean that the particle can enter the earth;
                   !   there is possibility so you have to examine.
                   !   if 0, there is no possibility so you may discard it.
!
!      Note: This system cannot treat multiple locations at a time.
!
!
      real*8 rigcut

      real*8  azm

      integer i1, i2, i3

      azm = azmin
      if(azm .lt. 0.) then
         azm = azm + 360.
      else
         azm = mod(azm, 360.d0)
      endif

      if(rdatafmt .eq. 1 ) then
         call k4ptdi(rigcuttbl, azmsize, zensize, azmsize, minazm, 
     *        minzen, dazm,  dzen,  azm, zen, rigcut)
         if(rig .lt. rigcut) then
            prob = 0.
         else
            prob = 1.
         endif
      elseif(rdatafmt .le. 4) then
!         i3 = log10(rig/MinRig)/LogDRig + 1
!         i3 <= log10(rig)*20 < i3+1 (i3=0,1..)
!         logrig = log10(rig)/LogDRig
!         i3 = logrig + 1 
!         dlogrig=(logrig - (i3-1))*LogDRig 

         i3 = log10(rig/minrig)/logdrig + 1
         if(i3 .lt. 1) then
            prob = 0.
            goto 100
         endif
         if(i3 .ge. rigsize) then
            prob = 1.
            goto 100
         endif
         if(rdatafmt .eq. 4) then
            i2=1
         else
            i2 =  azm/dazm + 1
            i2 =min(azmsize, i2)
         endif
         i1 = (zen-zenmax)/dzen +1     ! DZen < 0
         i1 = min(zensize, i1)
         call cgetrigprob(
     *    rigtbl2, zensize, azmsize, rigsize, i1, i2, i3,
     *    prob ) 
         prob = min(prob, 1.d0)
      elseif( rdatafmt .eq. 5) then
!           essentiall the same as 1 case. but rather should not
!           make interpolation
         call k2dtblv(rigcuttbl, azmsize, zensize, azmsize, minazm, 
     *        minzen, dazm,  dzen,  azm, zen, rigcut)
         if(rig .lt. rigcut) then
            prob = 0.
         else
            prob = 1.
         endif
      else
         call cerrormsg('format # error in rigidity table',0)
      endif
 100  continue

Here is the call graph for this function:

Here is the caller graph for this function:

crigcut0()

subroutine crigcut0

(

character*(*)

file

)

Definition at line 166 of file crigCut.f.

References cerrormsg(), copenf(), crigcutfmt1(), crigcutfmt2(), crigcutfmt4(), and crigcutfmt5().

Referenced by __drawrigcut.f__(), __drawrigcut2.f__(), and cbeginrun().

      implicit none
#include "Zmanagerp.h"
#include "ZrigCut.h"

      character*(*) file
      character*192 msg
      integer icon

      call copenf(tempdev, file, icon)
      if(icon .ne. 0) then
         write(msg, *) ' file specification error '
         call cerrormsg(msg, 0)
      endif
!
      read(tempdev, '(1x,i1)') rdatafmt   ! format 1 or 2 or 3 or 4
      if(rdatafmt .eq. 1) then
         call crigcutfmt1
      elseif(rdatafmt .eq. 2 .or. rdatafmt .eq. 3) then
         call crigcutfmt2
      elseif(rdatafmt .eq. 4) then
         call crigcutfmt4
      elseif(rdatafmt .eq. 5) then
         call crigcutfmt5
      else
         write(msg, *) 'rigidity cut data format =',rdatafmt,
     *   ' ivalid' 
         call cerrormsg(msg, 0)
      endif

Here is the call graph for this function:

Here is the caller graph for this function:

crigcutfmt1()

subroutine crigcutfmt1

(

)

Definition at line 197 of file crigCut.f.

References cerrormsg(), cfillrigcut(), creadrigcut(), cskipcomment(), and d0.

Referenced by crigcut0().

      implicit none
#include "Zmanagerp.h"
#include "ZrigCut.h"


      character*192 msg

      integer icon

      read(tempdev, *)   ! skip
      read(tempdev, *) place,  latit,  longi, magdec,  azmvalue,
     *            dazm, azmsize, zenvalue, dzen, zensize

!           skip until -------------- line
      call cskipcomment(tempdev, icon)
      if(icon .ne. 0) then
          call cerrormsg('cut-off table data has no ---- line', 0)
      endif
!          see if azimuthal angle range is from 0 to 360 or from 0 to 360-DAzm

      if(azmvalue .eq. 'deg') then
!               assume min of Azimuthal angle is  0 deg.
          minazm = 0.
          if( (azmsize - 1)* dazm  .lt. (360.d0- dazm* 0.1d0) ) then
!                 fill the last col.  by the first col.

             call creadrigcut(tempdev, rigcuttbl, azmsize, zensize, 
     *       azmsize+1)
             azmsize = azmsize + 1
             call cfillrigcut(rigcuttbl, azmsize, zensize)
          else
             call creadrigcut(tempdev, rigcuttbl, azmsize, zensize,
     *        azmsize)
          endif

       else

          write(msg, *)
     *    ' Azimuthal angle unit must be deg for rigidity cut table'
          call cerrormsg(msg, 1)
          write(msg, *) ' But it is ', azmvalue
          call cerrormsg(msg, 0)
       endif

       close(tempdev)

       if(zenvalue .eq. 'cos') then
!              DZen should be negative and min is 1.0
          minzen = 1.0
          if(dzen .ge. 0.) then
             write(msg, *)
     *        ' step of Zenith angle for rigidity cut should be < 0'
             call cerrormsg(msg, 1)
             write(msg, *) ' because you give it in cos value'
             call cerrormsg(msg, 0)
          endif
       endif
!
!                                   
       write(msg,*) 'Rigidity cut-off table has been read:',
     * ' place=',place,' latitute=',latit, ' longitude=',longi,
     * ' mag. dec=', magdec
       call cerrormsg(msg, 1)

Here is the call graph for this function:

Here is the caller graph for this function:

crigcutfmt2()

subroutine crigcutfmt2

(

)

Definition at line 308 of file crigCut.f.

References cerrormsg(), crdrigcut2(), and cskipcomment().

Referenced by crigcut0().

      implicit none
#include "Zmanagerp.h"
#include "ZrigCut.h"

      character*192 msg

      integer icon


      read(tempdev, *)   ! skip
      read(tempdev, *) 
     *  place,  latit,  longi, magdec, zenvalue, zenmax, dzen, zensize,
     *  azmvalue, minazm, dazm, azmsize, minrig, logdrig, rigsize


!           skip until -------------- line
      call cskipcomment(tempdev, icon)
      if(icon .ne. 0) then
          call cerrormsg('cut-off table data has no ---- line', 0)
      endif
      call crdrigcut2(
     *       rdatafmt, tempdev, rigtbl2, zensize, azmsize, rigsize)
      close(tempdev)

      if(zenvalue .eq. 'cos') then
!           DZen should be negative
         if(dzen .ge. 0.) then
            write(msg, *)
     *        ' step of Zenith angle for rigidity cut should be < 0'
            call cerrormsg(msg, 1)
            write(msg, *) ' because you give it in cos value'
            call cerrormsg(msg, 0)
         endif
      endif
!
!                                   
       write(msg,*) 'New rigidity cut-off table has been read:',
     * ' place=',place,' latitute=',latit, ' longitude=',longi,
     * ' mag. dec=', magdec
       call cerrormsg(msg, 1)

Here is the call graph for this function:

Here is the caller graph for this function:

crigcutfmt4()

subroutine crigcutfmt4

(

)

Definition at line 351 of file crigCut.f.

References cerrormsg(), crdrigcut4(), and cskipcomment().

Referenced by crigcut0().

      implicit none
#include "Zmanagerp.h"
#include "ZrigCut.h"

      character*192 msg

      integer icon

!      real pw

      read(tempdev, *)   ! skip
      read(tempdev, *) 
     *  place,  latit,  longi, magdec, zenvalue, zenmax, dzen, zensize,
     *  minrig, logdrig, rigsize

      azmsize=1


!           skip until -------------- line
      call cskipcomment(tempdev, icon)
      if(icon .ne. 0) then
          call cerrormsg('cut-off table data has no ---- line', 0)
      endif
      call crdrigcut4(
     *       rdatafmt, tempdev, minrig, rigtbl2, zensize, rigsize)
      close(tempdev)

      if(zenvalue .eq. 'cos') then
!           DZen should be negative
         if(dzen .ge. 0.) then
            write(msg, *)
     *        ' step of Zenith angle for rigidity cut should be < 0'
            call cerrormsg(msg, 1)
            write(msg, *) ' because you give it in cos value'
            call cerrormsg(msg, 0)
         endif
      endif
!
       write(msg,*) 'New rigidity cut-off table (fmt4) has been read:',
     * ' place=',place,' latitute=',latit, ' longitude=',longi,
     * ' mag. dec=', magdec
       call cerrormsg(msg, 1)

Here is the call graph for this function:

Here is the caller graph for this function:

crigcutfmt5()

subroutine crigcutfmt5

(

)

Definition at line 485 of file crigCut.f.

References cerrormsg(), cfillrigcut(), creadrigcut(), cskipcomment(), and d0.

Referenced by crigcut0().

      implicit none
#include "Zmanagerp.h"
#include "ZrigCut.h"

      character*192 msg

      integer icon


      read(tempdev, *)   ! skip
      read(tempdev, *)
     *  zenvalue, zenmax, dzen, zensize,
     *  azmvalue, minazm, dazm, azmsize


!           skip until -------------- line
      call cskipcomment(tempdev, icon)
      if(icon .ne. 0) then
          call cerrormsg('cut-off table data has no ---- line', 0)
      endif

      if(azmvalue .eq. 'deg') then
!               assume min of logitudinal angle is  0 deg.
          minazm = 0.
          if( (azmsize - 1)* dazm  .lt. (360.d0- dazm* 0.1d0) ) then
!                 fill the last col.  by the first col.
             call creadrigcut(tempdev, rigcuttbl, azmsize, zensize, 
     *       azmsize+1)
             azmsize = azmsize + 1
             call cfillrigcut(rigcuttbl, azmsize, zensize)
          else
             call creadrigcut(tempdev, rigcuttbl, azmsize, zensize,
     *        azmsize)
          endif
          close(tempdev)
       else
          write(msg, *)
     *    ' Azimuthal angle unit must be deg for rigidity cut table'
          call cerrormsg(msg, 1)
          write(msg, *) ' But it is ', azmvalue
          call cerrormsg(msg, 0)
       endif


       if(zenvalue .eq. 'cos') then
!              DZen should be negative and min is 1.0
          minzen = 1.0
          if(dzen .ge. 0.) then
             write(msg, *)
     *        ' step of Zenith angle for rigidity cut should be < 0'
             call cerrormsg(msg, 1)
             write(msg, *) ' because you give it in cos value'
             call cerrormsg(msg, 0)
          endif
       endif
!
!                                   
       write(msg,*) 'Rough rigidity cut-off table has been read:'
       call cerrormsg(msg, 1)

Here is the call graph for this function:

Here is the caller graph for this function:

k2dtblv()

subroutine k2dtblv	(	real*8, dimension(adj, ys)	tbl,
		integer	xs,
		integer	ys,
		integer	adj,
		real*8	xm,
		real*8	ym,
		real*8	dx,
		real*8	dy,
		real*8	x,
		real*8	y,
		real*8	ans
	)

Definition at line 97 of file crigCut.f.

Referenced by crigcut().

      integer xs
      integer ys
      integer adj
      real*8 xm
      real*8 ym
      real*8 dx
      real*8 dy
      real*8 x, y
      real*8 ans
      real*8 tbl(adj, ys)
      
      integer i, j
      
      i = (x-xm)/dx + 1
      j = (y-ym)/dy + 1
      ans = tbl(i,j)

Here is the caller graph for this function: