package scientific.stats.continuous

import std.math.*
import std.unittest.*
import std.unittest.testmacro.*

import scientific.numbers.*
import scientific.stats.random.*

public class Dstats {
    var mu: Float64
    var variance: Float64
    var skew: Float64
    var kurtosis: Float64

    init(mu: Float64, variance: Float64, skew: Float64, kurtosis: Float64) {
        this.mu = mu
        this.variance = variance
        this.skew = skew
        this.kurtosis = kurtosis
    }
}

/*
 * Log of Probability density function
 */
public func arcsineLogPDF(x: Float64,  loc!: Float64 = 0.0, scale!: Float64 = 1.0): Float64 {

    let y = (x - loc) / scale

    if(y <= 0.0 || y >= 1.0) {
        throw IllegalArgumentException("arcsineLogPDF: input value x out of bound.")
    }

    if(scale == 0.0) {
        throw IllegalArgumentException("arcsineLogPDF: scale cannot be 0.")
    }

    return - log(Float64.getPI()) - 0.5 * log(y) - 0.5 * log(1.0 - y) - log(scale)
}

/*
 * Probability density function
 */
public func arcsinePDF(x: Float64, loc!: Float64 = 0.0, scale!: Float64 = 1.0): Float64 {

    let y = (x - loc) / scale

    if(y <= 0.0 || y >= 1.0) {
        throw IllegalArgumentException("arcsinePDF: input value x out of bound.")
    }

    if(scale == 0.0) {
        throw IllegalArgumentException("arcsinePDF: scale cannot be 0.")
    }

    return exp(arcsineLogPDF(y, loc: loc, scale: scale))
}


/*
 * Cumulative probability density function
 */
public func arcsineCDF(x: Float64, loc!: Float64 = 0.0, scale!: Float64 = 1.0): Float64 {

    let y = (x - loc) / scale

    if(y <= 0.0 || y >= 1.0) {
        throw IllegalArgumentException("arcsineCDF: input value x out of bound.")
    }

    return 2.0 / Float64.getPI() * asin(sqrt(y))
}

/*
 * Log of cumulative probability density function
 */
public func arcsineLogCDF(x: Float64, loc!: Float64 = 0.0, scale!: Float64 = 1.0): Float64 {

    let y = (x - loc) / scale

    if(y <= 0.0 || y >= 1.0) {
        throw IllegalArgumentException("arcsineLogCDF: input value x out of bound.")
    }

    let temp = arcsineCDF(x, loc: loc, scale: scale)

    if(temp < 0.000001) {
        throw IllegalArgumentException("argusLogCDF: return-value too small.")
    }

    return log(temp)
}

/*
 * ppf
 */
public func arcsinePPF(q: Float64, loc!: Float64 = 0.0, scale!: Float64 = 1.0): Float64 {

    if(q <= 0.0 || q >= 1.0) {
        throw IllegalArgumentException("arcsinePPF: ppf point q out of bound.")
    }

    let temp = sin(Float64.getPI() * 0.5 * q)
    return temp * temp * scale + loc
}


/*
 * stats info: (mu, var, skew, kurtosis)
 */
public func arcsineStats(): Dstats {
    return Dstats(0.5, 1.0 / 8.0, 0.0, - 3.0 / 2.0)
}


/*
 * entropy of the RV
 */
public func arcsineEntropy(): Float64 {
    return -0.24156447527049044468
}


@Test
public class TestArcsine {
    @TestCase
    func testArcsine(): Unit {
        @Assert(approxEqual(arcsineLogPDF(1.2, loc: 1.0, scale: 2.0), -0.6339042620834093, atol:1e-13))
        @Assert(approxEqual(arcsineLogCDF(1.2, loc: 1.0, scale: 2.0), -1.5855614147303536, atol:1e-13))
        @Assert(approxEqual(arcsinePPF(0.2, loc: 1.0, scale: 2.0),     1.1909830056250525, atol:1e-13))
    }
}