diff --git a/NAMESPACE b/NAMESPACE
index b9fde9d..eeb0b53 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -15,6 +15,7 @@ export(dscaleGAL)
 export(mGAL)
 export(mean.variance.GMM.gradient.GAL)
 export(mgfGAL)
+export(pGAL)
 export(pnormapproxEsscherLM)
 export(psaddleapproxEsscherGAL)
 export(psaddleapproxGAL)
diff --git a/R/dGAL.r b/R/dGAL.r
index 664d537..5159627 100644
--- a/R/dGAL.r
+++ b/R/dGAL.r
@@ -1,3 +1,10 @@
+# Density function for the GAL distribution
+# 
+# Author: François Pelletier
+#
+# LGPL 3.0
+###############################################################################
+
 #' Density function for the GAL distribution
 #' @param x vector of quantiles
 #' @param param Parameter vector
diff --git a/R/pGAL.R b/R/pGAL.R
new file mode 100644
index 0000000..510108c
--- /dev/null
+++ b/R/pGAL.R
@@ -0,0 +1,21 @@
+# Distribution function for the GAL distribution
+# 
+# Author: François Pelletier
+#
+# LGPL 3.0
+###############################################################################
+
+#' Distribution function for the GAL distribution
+#' @param x vector of quantiles
+#' @param param Parameter vector
+#' @param type Choose between "mu" or "kappa" parametrization
+#' @param log Logical for log-parameters
+#' @return distribution at quantile x
+#' @export pGAL
+#' @author Francois Pelletier
+pGAL <- function(x,param,type="mu",log=FALSE)
+{
+	x.l <- as.list(x)
+	integrate1 <- function(x,param,type,log) integrate(dGAL,-Inf,x,param=param,type=type,log=log)$value
+	unlist(lapply(x.l,integrate1,param,type,log))
+}
diff --git a/R/psaddleapproxGAL.R b/R/psaddleapproxGAL.R
index 5d27351..f5fd036 100644
--- a/R/psaddleapproxGAL.R
+++ b/R/psaddleapproxGAL.R
@@ -22,8 +22,8 @@ psaddleapproxGAL <- function(x,param,eval.time=1,type="mu",log=FALSE)
 	u <- s * sqrt(diffcgfGAL(s,2,param,eval.time,type,log))
 	w <- sign(s)*sqrt(2*(s*x-cgfGAL(s,param,type,log)))
 	
-	(x==mGAL(param,1,type,log)) * (1/2 + diffcgfGAL(0,3,param,eval.time,type,log)/
+	(x==round(mGAL(param,1,type,log),4)) * (1/2 + diffcgfGAL(0,3,param,eval.time,type,log)/
 				(6*sqrt(2*pi)*diffcgfGAL(0,2,param,eval.time,type,log)^(3/2))) + 
-			(x!=mGAL(param,1,type,log)) * pnorm(w)+dnorm(w)*(1/w-1/u)
+			(x!=round(mGAL(param,1,type,log),4)) * (pnorm(w)+dnorm(w)*(1/w-1/u))
 			
 }
diff --git a/man/pGAL.Rd b/man/pGAL.Rd
new file mode 100644
index 0000000..3f5559b
--- /dev/null
+++ b/man/pGAL.Rd
@@ -0,0 +1,26 @@
+\name{pGAL}
+\alias{pGAL}
+\title{Distribution function for the GAL distribution}
+\usage{
+pGAL(x, param, type = "mu", log = FALSE)
+}
+\arguments{
+  \item{x}{vector of quantiles}
+
+  \item{param}{Parameter vector}
+
+  \item{type}{Choose between "mu" or "kappa"
+  parametrization}
+
+  \item{log}{Logical for log-parameters}
+}
+\value{
+distribution at quantile x
+}
+\description{
+Distribution function for the GAL distribution
+}
+\author{
+Francois Pelletier
+}
+