cubehelix_colormap Subroutine

public pure subroutine cubehelix_colormap(map, nlev, varargs)

Based on the public domain FORTRAN 77 subroutine published by D.A. Green: Green, D. A., 2011, Bulletin of the Astronomical Society of India, Vol.39, p.289 For more information on the parameters of cubehelix, see his page: https://www.mrao.cam.ac.uk/~dag/CUBEHELIX/

Arguments

Type IntentOptional Attributes Name
integer, intent(out), dimension(:,:), allocatable :: map
integer, intent(in) :: nlev
real(kind=wp), intent(in), optional, dimension(:) :: varargs

Calls

proc~~cubehelix_colormap~~CallsGraph proc~cubehelix_colormap cubehelix_colormap local local proc~cubehelix_colormap->local

Called by

proc~~cubehelix_colormap~~CalledByGraph proc~cubehelix_colormap cubehelix_colormap proc~set Colormap%set proc~set->proc~cubehelix_colormap proc~test_038 test_038 proc~test_038->proc~set proc~test_039 test_039 proc~test_039->proc~set proc~test_040 test_040 proc~test_040->proc~set proc~test_045 test_045 proc~test_045->proc~set proc~test_048 test_048 proc~test_048->proc~set proc~test_056 test_056 proc~test_056->proc~set proc~test_057 test_057 proc~test_057->proc~set proc~test_058 test_058 proc~test_058->proc~set proc~test_059 test_059 proc~test_059->proc~set proc~test_061 test_061 proc~test_061->proc~set proc~test_062 test_062 proc~test_062->proc~set proc~test_063 test_063 proc~test_063->proc~set proc~test_068 test_068 proc~test_068->proc~set proc~test_074 test_074 proc~test_074->proc~set proc~test_076 test_076 proc~test_076->proc~set proc~test_077 test_077 proc~test_077->proc~set proc~test_078 test_078 proc~test_078->proc~set proc~test_079 test_079 proc~test_079->proc~set proc~test_080 test_080 proc~test_080->proc~set proc~test_082 test_082 proc~test_082->proc~set proc~test_083 test_083 proc~test_083->proc~set proc~test_085 test_085 proc~test_085->proc~set program~demo demo program~demo->proc~set program~demo_reverse demo_reverse program~demo_reverse->proc~set program~extract extract program~extract->proc~set program~modify modify program~modify->proc~set program~check check program~check->proc~test_038 program~check->proc~test_039 program~check->proc~test_040 program~check->proc~test_045 program~check->proc~test_048 program~check->proc~test_056 program~check->proc~test_057 program~check->proc~test_058 program~check->proc~test_059 program~check->proc~test_061 program~check->proc~test_062 program~check->proc~test_063 program~check->proc~test_068 program~check->proc~test_074 program~check->proc~test_076 program~check->proc~test_077 program~check->proc~test_078 program~check->proc~test_079 program~check->proc~test_080 program~check->proc~test_082 program~check->proc~test_083 program~check->proc~test_085

Source Code

    pure subroutine cubehelix_colormap(map, nlev, varargs)
        integer, dimension(:,:), allocatable, intent(out) :: map
        integer, intent(in) :: nlev
        real(wp), dimension(:), intent(in), optional :: varargs
        integer  :: i
        real(wp) :: start, rots, hue, gamma
        real(wp) :: fract, angle, amp

        if (nlev < 2) error stop "ERROR: cubehelix nlev must be >= 2"

        if (present(varargs)) then
            if (size(varargs) /= 4) error stop "ERROR: cubehelix varargs(:) must have 4 values"
            start = varargs(1)
            rots  = varargs(2)
            hue   = varargs(3)
            gamma = varargs(4)
        else
            ! Default values:
            start = 0.5_wp
            rots  = -1.5_wp
            hue   = 1.0_wp
            gamma = 1.0_wp
        end if

        allocate(map(0:nlev-1, 1:3))

        do concurrent (i = 0:nlev-1) local(fract, angle, amp)
            fract = real(i, kind=wp) / (nlev-1)
            angle = 2*pi * (start/3 + 1 + rots*fract)
            fract = fract ** gamma
            amp   = hue * fract * (1-fract)/2

            map(i, 1) = nint(255*(fract + amp*(-0.14861_wp*cos(angle) + 1.78277_wp*sin(angle))))
            map(i, 2) = nint(255*(fract + amp*(-0.29227_wp*cos(angle) - 0.90649_wp*sin(angle))))
            map(i, 3) = nint(255*(fract + amp*(+1.97294_wp*cos(angle))))
        end do
    end subroutine cubehelix_colormap