init_spectral Function

    type(spectral) function init_spectral(mus, mua, hgg, n, flux) result(res)

        real(kind=wp), allocatable, intent(in) :: mus(:, :), mua(:, :), hgg(:, :), n(:, :), flux(:, :)
        
        real(kind=wp) :: wave, tmp

        !setup cdfs
        res%flux = piecewise1D(flux)

        res%mus_a = piecewise1D(mus)
        res%mua_a = piecewise1D(mua)
        res%hgg_a = piecewise1D(hgg)
        res%n_a = piecewise1D(n)

        !sample wavelength so we can sample from other optical properties at the correct points
        call res%flux%sample(wave,  tmp)
        ! sample optical properties
        call res%mus_a%sample(res%mus, wave)
        call res%mua_a%sample(res%mua, wave)
        call res%hgg_a%sample(res%hgg, wave)
        res%g2 = res%hgg**2
        call res%n_a%sample(res%n, wave)

        res%kappa =  res%mus + res%mua
        if(res%mua < 1e-9_wp)then          
            res%albedo = 1.
        else
            res%albedo = res%mus / res%kappa
        end if
    end function init_spectral