Module helixfit_cell_mod

module helixfit_cell_mod

        ! Uses
    use precision_mod
    use namelist_mod
    use physicalconstants_mod
    use track_mod
    use helixfit_common_mod
    use helixfit_hits_mod

        ! Variables
    logical, private, SAVE, dimension (kMaxPlanes) :: haveTable6 = .false.
    logical, private, SAVE, dimension (kMaxPlanes) :: haveTable7 = .false.
    logical, private, SAVE, dimension (kMaxPlanes) :: haveTable8 = .false.

        ! Subroutines and functions
    public subroutine ResetCellTables ()
    private subroutine ClosestToWire (V, chit, X)
    public function DistanceToWire (V, chit) result (dist)
    private function getDriftDistance (v, p, dir, chit, debug) result (returnvalue)
    public function gotoHit (hits, iplane, ihit, V1, V2, V3, V, P, Q, K, enableSTR, minpos, mindrift, mintime, tmean, trms, isInside, wirepos) result (ierror)
    public subroutine testGotoHit (hits)
    public function getHitDrift (hits, iplane, ihit, trkpos, trkmom, drifttime, hitpos3) result (driftdist)
    private subroutine testIsochrone (r, phi)
    private subroutine testGetCellMinTimeOne (iplane, r, xphi, minpos, mintime, debug)
    public subroutine testGetCellMinTime ()
    private subroutine testGetCellMinTimeCorrection ()
    private subroutine testTableDTT ()
    private subroutine testTable6 ()
    private subroutine getCellCoord (planedir, wirepos, uvzpos, cellpos)
    private function getCellTime (iplane, pos) result (time)
    private function getCellTime1 (iplane, pos) result (time)
    private subroutine getCellMinDistQI3 (pos0, pos1, pos2, pos3, tr, minpos, debug)
    private subroutine getCellMinTimeGoldenSection (iplane, tr, posEnter, posExit, precz, prect, pos1, time1, pos2, time2, debug)
    private subroutine getCellMinTimeGrid (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    private subroutine getCellMinTime (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    private subroutine getIsochroneTable (iplane, t, phi, r, drdphi)
    private subroutine getIsochroneTangentTable (iplane, phi, r, drdphi, duvdz)
    private subroutine getCellMinTimeIsochrone (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    private subroutine getCellMinTimeLinear (iplane, r, xphi, minpos, mintime, debug)
    private subroutine getCellMinTimeCorrection (iplane, r, phi, cellz, debug)
    private subroutine getCellMinTimeCorrectionTable (iplane, r, phi, cellz, debug)
    private subroutine getQuad (pos1, pos2, pos3, tr)
    private function zpos (tr, z) result (pos)
    private function ztangent (tr, z)
    private function uvpos (tr, zref, uv) result (pos)
    private function fitQuadratic (s1, s2, s3, t1, t2, t3, c0, c1, c2) result (retvalue)
    private function myReal8 (a4) result (a8)
    private function jsign (a) result (ival)
    private function V3mag (a)
    private function sqr (a) result (b)
    private function clip (x, c) result (y)
    private function prec (x, p) result (y)

end module helixfit_cell_mod

Description of Variables

logical, private, SAVE, dimension (kMaxPlanes) :: haveTable6 = .false.
 isochrone table

logical, private, SAVE, dimension (kMaxPlanes) :: haveTable7 = .false.
 isochrone table

logical, private, SAVE, dimension (kMaxPlanes) :: haveTable8 = .false.
 isochrone table

Description of Subroutines and Functions

public subroutine ResetCellTables ()
end subroutine ResetCellTables

private subroutine ClosestToWire (V, chit, X)
    real (kind=R8), INTENT(in), dimension (3) :: V
            
 Description: calculate and return point closest to wire

 Arguments:

 point
    type (hit_type), INTENT(in) :: chit
             hit data
    real (kind=R8), INTENT(out), dimension (3) :: X
             point closest to the wire
end subroutine ClosestToWire

public function DistanceToWire (V, chit) result (dist)
    real (kind=R8), INTENT(in), dimension (3) :: V
            
 Arguments:

 point coordinates
    type (hit_type), INTENT(in) :: chit
             wire data
    real (kind=R8) :: dist
            
 Description: calculate distance from point to wire

 return value
    ! Calls: ClosestToWire
end function DistanceToWire

private function getDriftDistance (v, p, dir, chit, debug) result (returnvalue)
    real (kind=R8), INTENT(in), dimension (3) :: v
             track position
    real (kind=R8), INTENT(in), dimension (3) :: p
             track momentum
    real (kind=R8), INTENT(in), dimension (3) :: dir
    type (hit_type), INTENT(in) :: chit
             hit data
    logical, INTENT(in) :: debug
    real (kind=R8) :: returnvalue
             return value
    ! Calls: closestToWire
end function getDriftDistance

public function gotoHit (hits, iplane, ihit, V1, V2, V3, V, P, Q, K, enableSTR, minpos, mindrift, mintime, tmean, trms, isInside, wirepos) result (ierror)
    type (hits_type), INTENT(in) :: hits
            
 Description: find predicted position of this hit

 Arguments:

 hits
    integer, INTENT(in) :: iplane
    integer, INTENT(in) :: ihit
    real (kind=R8), INTENT(in), dimension (3) :: V1
             track position at plane entrance
    real (kind=R8), INTENT(in), dimension (3) :: V2
             track position at plane z
    real (kind=R8), INTENT(in), dimension (3) :: V3
             track position at plane exit
    real (kind=R8), INTENT(inout), dimension (3) :: V
             new track position
    real (kind=R8), INTENT(inout), dimension (3) :: P
             new track position
    integer, INTENT(in) :: Q
             track charge
    real (kind=R8), INTENT(in) :: K
             track energy loss
    logical, INTENT(in) :: enableSTR
             enable use of drift tables
    real (kind=R8), INTENT(out), dimension (2) :: minpos
             cell position
    real (kind=R8), INTENT(out) :: mindrift
             signed (left-right) drift distance
    real (kind=R8), INTENT(out) :: mintime
             unsigned drift time
    real (kind=R8), INTENT(out) :: tmean
             clustered drift time mean
    real (kind=R8), INTENT(out) :: trms
             clustered drift time rms
    logical, INTENT(out) :: isInside
             .true. is hit is inside the drift cell
    real (kind=R8), INTENT(out), dimension (3) :: wirepos
             hit wire position
    integer :: ierror
             return value
    ! Calls: getCellCoord, getCellMinDistQI3, getCellMinTime, getCellSizeZ, getQuad
end function gotoHit

public subroutine testGotoHit (hits)
    type (hits_type), INTENT(in) :: hits
            
 Description: test the gotoHit() function

 Arguments:

 hits
end subroutine testGotoHit

public function getHitDrift (hits, iplane, ihit, trkpos, trkmom, drifttime, hitpos3) result (driftdist)
    type (hits_type), INTENT(in) :: hits
            
 Description: using the trkpos hit position
    and the measured drift time, calculate
    the measured drift distance.

 hits record
    integer, INTENT(in) :: iplane
    integer, INTENT(in) :: ihit
    real (kind=R8), INTENT(in), dimension (3) :: trkpos
    real (kind=R8), INTENT(in), dimension (3) :: trkmom
    real (kind=R8), INTENT(in) :: drifttime
             drift time
    real (kind=R8), INTENT(out), dimension (3) :: hitpos3
             hit position
    real (kind=R8) :: driftdist
             return value
    ! Calls: getCellCoord, getIsochroneTable
end function getHitDrift

private subroutine testIsochrone (r, phi)
    real (kind=R8), intent(IN) :: r
    real (kind=R8), intent(IN) :: phi
    ! Calls: getIsochroneTable
end subroutine testIsochrone

private subroutine testGetCellMinTimeOne (iplane, r, xphi, minpos, mintime, debug)
    integer, intent(IN) :: iplane
    real (kind=R8), intent(IN) :: r
    real (kind=R8), intent(IN) :: xphi
    real (kind=R8), intent(OUT), dimension (2) :: minpos
    real (kind=R8), intent(OUT) :: mintime
    logical, intent(IN) :: debug
    ! Calls: getCellMinTime, getCellSizeZ, getQuad
end subroutine testGetCellMinTimeOne

public subroutine testGetCellMinTime ()
    ! Calls: testGetCellMinTimeCorrection, testGetCellMinTimeOne
end subroutine testGetCellMinTime

private subroutine testGetCellMinTimeCorrection ()
    ! Calls: getCellMinTimeCorrectionTable
end subroutine testGetCellMinTimeCorrection

private subroutine testTableDTT ()
end subroutine testTableDTT

private subroutine testTable6 ()
    ! Calls: getIsochroneTable
end subroutine testTable6

private subroutine getCellCoord (planedir, wirepos, uvzpos, cellpos)
    real (kind=R8), INTENT(in), dimension (3) :: planedir
    real (kind=R8), INTENT(in), dimension (3) :: wirepos
    real (kind=R8), INTENT(in), dimension (3) :: uvzpos
    real (kind=R8), INTENT(out), dimension (2) :: cellpos
end subroutine getCellCoord

private function getCellTime (iplane, pos) result (time)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (2) :: pos
    real (kind=R8) :: time
end function getCellTime
 call normal interpolated drift tables (quadratic interpolation)

private function getCellTime1 (iplane, pos) result (time)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (2) :: pos
    real (kind=R8) :: time
end function getCellTime1
 call linearly interpolated drift tables

private subroutine getCellMinDistQI3 (pos0, pos1, pos2, pos3, tr, minpos, debug)
    real (kind=R8), INTENT(in), dimension (2) :: pos0
            
 Description: on a trajectory passing through pos1, pos2 and pos3
 find the point closest to pos0 and return it's coordinates
 in minpos.
    real (kind=R8), INTENT(in), dimension (2) :: pos1
    real (kind=R8), INTENT(in), dimension (2) :: pos2
    real (kind=R8), INTENT(in), dimension (2) :: pos3
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(out), dimension (2) :: minpos
    logical, INTENT(in) :: debug
end subroutine getCellMinDistQI3

private subroutine getCellMinTimeGoldenSection (iplane, tr, posEnter, posExit, precz, prect, pos1, time1, pos2, time2, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in), dimension (2) :: posEnter
    real (kind=R8), INTENT(in), dimension (2) :: posExit
    real (kind=R8), INTENT(in) :: precz
    real (kind=R8), INTENT(in) :: prect
    real (kind=R8), INTENT(out), dimension (2) :: pos1
    real (kind=R8), INTENT(out) :: time1
    real (kind=R8), INTENT(out), dimension (2) :: pos2
    real (kind=R8), INTENT(out) :: time2
    logical, INTENT(in) :: debug
end subroutine getCellMinTimeGoldenSection

private subroutine getCellMinTimeGrid (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in), dimension (2) :: posEnter
    real (kind=R8), INTENT(in), dimension (2) :: posExit
    real (kind=R8), INTENT(out), dimension (2) :: minpos
    real (kind=R8), INTENT(out) :: mintime
    logical, INTENT(in) :: debug
end subroutine getCellMinTimeGrid

private subroutine getCellMinTime (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in), dimension (2) :: posEnter
    real (kind=R8), INTENT(in), dimension (2) :: posExit
    real (kind=R8), INTENT(out), dimension (2) :: minpos
    real (kind=R8), INTENT(out) :: mintime
    logical, INTENT(in) :: debug
    ! Calls: getCellMinDistQI3, getCellMinTimeCorrectionTable
end subroutine getCellMinTime

private subroutine getIsochroneTable (iplane, t, phi, r, drdphi)
    integer, INTENT(in) :: iplane
             plane number
    real (kind=R8), INTENT(in) :: t
             find isochorone with time t, ns
    real (kind=R8), INTENT(in) :: phi
             angular direction, deg
    real (kind=R8), INTENT(out) :: r
             radius of the isochrone, cm
    real (kind=R8), INTENT(out) :: drdphi
             dr/dphi tangent to the isochrone
    ! Calls: buildTable6
end subroutine getIsochroneTable

private subroutine getIsochroneTangentTable (iplane, phi, r, drdphi, duvdz)
    integer, INTENT(in) :: iplane
             plane number
    real (kind=R8), INTENT(in) :: phi
             angular direction, deg
    real (kind=R8), INTENT(in) :: r
             radius of the isochrone, cm
    real (kind=R8), INTENT(out) :: drdphi
             dr/dphi tangent to the isochrone
    real (kind=R8), INTENT(out) :: duvdz
             duv/dz  tangent to the isochrone
    ! Calls: buildTable7
end subroutine getIsochroneTangentTable

private subroutine getCellMinTimeIsochrone (iplane, tr, posEnter, posExit, minpos, mintime, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in), dimension (2) :: posEnter
    real (kind=R8), INTENT(in), dimension (2) :: posExit
    real (kind=R8), INTENT(out), dimension (2) :: minpos
    real (kind=R8), INTENT(out) :: mintime
    logical, INTENT(in) :: debug
    ! Calls: bracketsSearch, die, f, fitIso2, getCellMinTimeGoldenSection, gridSearch, iso2search, linearSearch, probe, probeIso2, singularitySearch
end subroutine getCellMinTimeIsochrone

private subroutine getCellMinTimeLinear (iplane, r, xphi, minpos, mintime, debug)
    integer, intent(IN) :: iplane
    real (kind=R8), intent(IN) :: r
    real (kind=R8), intent(IN) :: xphi
    real (kind=R8), intent(OUT), dimension (2) :: minpos
    real (kind=R8), intent(OUT) :: mintime
    logical, intent(IN) :: debug
    ! Calls: getCellMinTimeIsochrone, getCellSizeZ, getQuad
end subroutine getCellMinTimeLinear

private subroutine getCellMinTimeCorrection (iplane, r, phi, cellz, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in) :: r
             impact parameter, cm
    real (kind=R8), INTENT(in) :: phi
             polar angle, deg
    real (kind=R8), INTENT(out) :: cellz
             z coordinate of point of minimum time
    logical, INTENT(in) :: debug
    ! Calls: getCellMinTimeLinear, getCellSizeZ
end subroutine getCellMinTimeCorrection

private subroutine getCellMinTimeCorrectionTable (iplane, r, phi, cellz, debug)
    integer, INTENT(in) :: iplane
    real (kind=R8), INTENT(in) :: r
             impact parameter, cm
    real (kind=R8), INTENT(in) :: phi
             polar angle, deg
    real (kind=R8), INTENT(out) :: cellz
             z coordinate of point of minimum time
    logical, INTENT(in) :: debug
    ! Calls: buildTable8
end subroutine getCellMinTimeCorrectionTable

private subroutine getQuad (pos1, pos2, pos3, tr)
    real (kind=R8), intent(IN), dimension (2) :: pos1
    real (kind=R8), intent(IN), dimension (2) :: pos2
    real (kind=R8), intent(IN), dimension (2) :: pos3
    real (kind=R8), intent(OUT), dimension (3) :: tr
    ! Calls: die
end subroutine getQuad
 approximate helical trajectory with a parabola

private function zpos (tr, z) result (pos)
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in) :: z
    real (kind=R8), dimension (2) :: pos
             return value
end function zpos
 position along the path

private function ztangent (tr, z)
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in) :: z
    real (kind=R8) :: ztangent
             return value
end function ztangent
 tangent to the path

private function uvpos (tr, zref, uv) result (pos)
    real (kind=R8), INTENT(in), dimension (3) :: tr
    real (kind=R8), INTENT(in) :: zref
    real (kind=R8), INTENT(in) :: uv
    real (kind=R8), dimension (2) :: pos
             return value
end function uvpos
 position along the path

private function fitQuadratic (s1, s2, s3, t1, t2, t3, c0, c1, c2) result (retvalue)
    real (kind=R8), INTENT(in) :: s1
    real (kind=R8), INTENT(in) :: s2
    real (kind=R8), INTENT(in) :: s3
    real (kind=R8), INTENT(in) :: t1
    real (kind=R8), INTENT(in) :: t2
    real (kind=R8), INTENT(in) :: t3
    real (kind=R8), INTENT(out) :: c0
    real (kind=R8), INTENT(out) :: c1
    real (kind=R8), INTENT(out) :: c2
    logical :: retvalue
             return value
end function fitQuadratic
 helper solver for quadratic equations

private function myReal8 (a4) result (a8)
    real, INTENT(in) :: a4
    real (kind=R8) :: a8
end function myReal8

private function jsign (a) result (ival)
    real (kind=R8), INTENT(in) :: a
            
 Description: return sign of a real(r8) number
 Returns: +1 if positive
          -1 if negative
           0 if zero
    integer :: ival
end function jsign

private function V3mag (a)
    real (kind=R8), INTENT(in), dimension (3) :: a
            
 Description: compute vector length
    real (kind=R8) :: v3mag
end function V3mag

private function sqr (a) result (b)
    real (kind=R8), INTENT(in) :: a
            
 Description: return square of a real(r8) number
    real (kind=R8) :: b
end function sqr

private function clip (x, c) result (y)
    real (kind=R8), INTENT(in) :: x
            
 Description: make abs(x) <= c
    real, INTENT(in) :: c
    real (kind=R8) :: y
             return value
end function clip

private function prec (x, p) result (y)
    real (kind=R8), intent(IN) :: x
             limit value of x to precision p
    real, intent(IN) :: p
    real (kind=R8) :: y
end function prec