008_Fisher_calculator_Cls_windows.f90

Go to the documentation of this file.

!----------------------------------------------------------------------------------------
!
! This file is part of CosmicFish.
!
! Copyright (C) 2015-2016 by the CosmicFish authors
!
! The CosmicFish code is free software;
! You can use it, redistribute it, and/or modify it under the terms
! of the GNU General Public License as published by the Free Software Foundation;
! either version 3 of the License, or (at your option) any later version.
! The full text of the license can be found in the file LICENSE at
! the top level of the CosmicFish distribution.
!
!----------------------------------------------------------------------------------------


!----------------------------------------------------------------------------------------


module fisher_calculator_cls_windows

    use precision
    use constants, only: const_pi
    use redshiftspacedata
    use cosmicfish_types

    implicit none

    type(cosmicfish_params), save :: window_fp

    private

    public init_camb_sources_windows, gaussian_window, bin_window, flat_window, flat_smooth_window

contains

    ! ---------------------------------------------------------------------------------------------
    subroutine init_camb_sources_windows( FP )

        implicit none

        Type(cosmicfish_params), intent(in)  :: FP

        window_fp = fp

        if ( fp%fisher_cls%window_type == 1 ) then
            window_function_calc => default_window
        else if ( fp%fisher_cls%window_type == 2 ) then
            window_function_calc => bin_window
        else if ( fp%fisher_cls%window_type == 3 ) then
            window_function_calc => flat_window
        else if ( fp%fisher_cls%window_type == 4 ) then
            window_function_calc => flat_smooth_window
        else
            window_function_calc => null()
        end if

    end subroutine init_camb_sources_windows

    ! ---------------------------------------------------------------------------------------------
    function gaussian_window(Win, z, winamp)

        implicit none

        Type(tredwin)        :: Win
        real(dl), intent(in) :: z
        real(dl)             :: winamp

        real(dl) :: gaussian_window

        real(dl) :: dz

        gaussian_window = 0._dl

        dz     = z - win%Redshift
        winamp =  exp(-0.5_dl*(dz/win%sigma)**2)
        gaussian_window = winamp/win%sigma/sqrt(2._dl*const_pi)

    end function gaussian_window

    ! ---------------------------------------------------------------------------------------------
    function bin_window(Win, z, winamp)

        implicit none

        Type(tredwin)        :: Win
        real(dl), intent(in) :: z
        real(dl)             :: winamp

        real(dl) :: bin_window

        real(dl) :: alpha, beta, z0, sigma, total_number_z

        alpha = window_fp%fisher_cls%window_alpha
        beta  = window_fp%fisher_cls%window_beta
        z0    = window_fp%fisher_cls%redshift_zero
        sigma = window_fp%fisher_cls%photoz_error*(1._dl+z)

        total_number_z = beta*z**alpha*exp(-(z/z0)**beta)/(z0**(alpha+1._dl))/gamma((alpha+1._dl)/beta)

        winamp = 0.5_dl*total_number_z*( &
            & +erf( (win%Redshift + win%sigma - z)/(sqrt(2._dl))/sigma ) &
            & -erf( (win%Redshift - z)/(sqrt(2._dl))/sigma ) )

        if ( winamp<1.d-10 ) winamp = 0._dl

        bin_window = winamp

    end function bin_window

    ! ---------------------------------------------------------------------------------------------
    function flat_window(Win, z, winamp)

        implicit none

        Type(tredwin)        :: Win
        real(dl), intent(in) :: z
        real(dl)             :: winamp

        real(dl) :: flat_window

        real(dl) :: alpha, beta, z0, sigma, total_number_z

        winamp = 0._dl

        if ( win%Redshift<z .and. z<win%Redshift + win%sigma ) then
            winamp = 1._dl
        end if

        flat_window = winamp

    end function flat_window

    ! ---------------------------------------------------------------------------------------------
    function flat_smooth_window(Win, z, winamp)

        implicit none

        Type(tredwin)        :: Win
        real(dl), intent(in) :: z
        real(dl)             :: winamp

        real(dl) :: flat_smooth_window

        real(dl) :: alpha, beta, z0, sigma, total_number_z

        winamp = 0._dl

        if ( win%Redshift<z .and. z<win%Redshift + win%sigma ) then
            winamp = 1._dl
        end if

        sigma = window_fp%fisher_cls%photoz_error*(1._dl+z)

        winamp = 0.5_dl*winamp*( &
            & +erf( (win%Redshift + win%sigma - z)/(sqrt(2._dl))/sigma ) &
            & -erf( (win%Redshift - z)/(sqrt(2._dl))/sigma ) )

        if ( winamp<1.d-10 ) winamp = 0._dl

        flat_smooth_window = winamp

    end function flat_smooth_window

end module fisher_calculator_cls_windows