# Roger Palay copyright 2016-02-15
# Saline, MI 48176
#
papfelton <- function( v, lower.tail=TRUE)
  { # compute the cumulative probability 
    # for the Apfelton distribution for vaules between 
    # -2 and 5 inclusive.
    # By default return P(X<v), but if lower.tail=FALSE
    # then return the upper tail
  
    # handle the extreme cases first
    if( v < -2 )
      {tail<-0}
    else if (v > 5 )
      {tail<-1}
    else
      {
       tail <- (v^5/5 - v^4 -5*v^3+29*v^2+128*v+122.4)/862.4
      }
    if( lower.tail)
       {return(tail)}
    else
       { return(1-tail)}
  }
  
qapfelton <- function( q, lower.tail=TRUE)
  { # inverse function of papfelton()
  if( q < 0 )
     { q <- 0 }
  else if ( q > 1 )
     {q <- 1 }
  if( !lower.tail )
     { q <- 1-q}
  # handle extreme cases
  if ( q == 0 ) {return(-2)}
  if ( q == 1 ) { return(5)}
  low <- -2
  high <- 5
  error <- 1
  while (error > 0.0000000001)
    { 
     next_v <- (low+high)/2
     next_p <- papfelton( next_v )
     if ( next_p <= q )
       { 
        low <- next_v
        error <- q - next_p
       }
     if ( next_p >= q )
     { 
       high <- next_v
       error <-  next_p - q
     }     
    }
  return( next_v )
  }