.packageName <- "NORMT3"
"dnormt3" <-
function(x, mean=0, sd=1){

if (sd==0)	{
	f <- rep(0, length(x))
	f[x==0] <- Inf
	}
else	{
	f <- sqrt(2)*normt3ip(mu=(x-mean)*sqrt(2)/sd, sigma=1)/sd
	sv <- abs(x - mean)*sqrt(2)/sd >= 37
	f[sv] <- 0
	}
f
}
"dst" <-
function(x,nu=3){
gamma( (nu+1)/2)/(gamma(nu/2)*sqrt(pi)*sqrt(nu-2)*( 1+ (x^2)/(nu-2))^((nu+1)/2))
}
"erf" <-
function(z) 1 - erfc(z)
"erfc" <-
function(z){

ans <- .C("IPerfcvec",
		x=as.double(Re(z)),
		y=as.double(Im(z)),
		ansx=as.double(Re(z)),
		ansy=as.double(Im(z)),
		n = as.integer(length(z)),
		error = as.integer(0),
		PACKAGE="NORMT3")
if (ans$error != 0)
	stop(paste("Error code from TOMS 680 was ", ans$error))

return(complex(real=ans$ansx, im=ans$ansy))

}
.First.lib <- function(lib,pkg)
{
   library.dynam("NORMT3",pkg,lib)
   cat("NORMT3 loaded\n")
}
"ic1" <-
function (p,d) 
{
cc <- complex(re=0, im=p/2)

(sqrt(pi)/2)*exp(-(p^2)/4)*(1-Re(erf(cc) - erf(cc-d)))

}
"is1" <-
function (p,d) 
{

cc <- complex(re=0, im=p/2)

(sqrt(pi)/2)*exp(-(p^2)/4)*Im(erf(cc-d))
}
"normt3ip" <-
function (mu, sigma) 
{
a <- 1/(sigma^2)
b <- sqrt(2)/sigma
d <- a/b
p <- mu*b

T1 <- ((1-a)*cos(mu*a) + p*b*sin(mu*a)/2)* ic1(p,d)

T2 <- ((1-a)*sin(mu*a) - p*b*cos(mu*a)/2)*is1(p,d)

T3 <- b*exp( - d^2)/2

total <- b*exp(a/2)*(T1+T2+T3)/pi
total
}
