-
Notifications
You must be signed in to change notification settings - Fork 103
/
YilmIndexVector.f95
65 lines (53 loc) · 1.67 KB
/
YilmIndexVector.f95
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
function YilmIndexVector(i, l, m)
!-------------------------------------------------------------------------------
!
! This function will give the index in a 1-dimensional array of the spherical
! harmonic coefficient corresponding to the element Cilm. The elements of the
! 1D vector array are packed according to (where positive m corresponds
! to i=1 (the cosine coefficients), and negative m corresponds to i=1
! (the sine coefficients))
!
! 0,0
! 1, 0; 1, 1; 1, -1
! 2,0 ; 2, 1; 2, 2; 2, -1; 2, -2
!
! This mapping is given by the function:
!
! YilmIndex = l**2 + (i-1)*l + m + 1
!
! Copyright (c) 2005-2019, SHTOOLS
! All rights reserved.
!
!-------------------------------------------------------------------------------
use ftypes
implicit none
integer(int32) :: YilmIndexVector
integer(int32), intent(in) :: i, l, m
if (i /= 1 .and. i /= 2) then
print*, "Error --- YilmIndexVector"
print*, "I must be 1 (for cosine terms) or 2 (for sine terms)."
print*, "I = ", i
stop
end if
if (l < 0) then
print*, "Error --- YilmIndexVector"
print*, "L must be positive."
print*, "L = ", l
stop
end if
if (m < 0 .or. m > l) then
print*, "Error --- YilmIndexVector"
print*, "M must be positive and less than L."
print*, "M = ", m
print*, "L = ", l
stop
end if
if (m == 0 .and. i == 2) then
print*, "Error --- YilmIndexVector"
print*, "When M = 0, I must be 1."
print*, "I = ", i
print*, "M = ", m
stop
end if
yilmindexvector = l**2 + (i-1)*l + m + 1
end function YilmIndexVector