source: LMDZ6/branches/DYNAMICO-conv-GC/tools/Max_diff_nc_with_lib/Jumble/Numerical/spherical.f90

Last change on this file was 1907, checked in by lguez, 11 years ago

Added a copyright property to every file of the distribution, except
for the fcm files (which have their own copyright). Use svn propget on
a file to see the copyright. For instance:

$ svn propget copyright libf/phylmd/physiq.F90
Name of program: LMDZ
Creation date: 1984
Version: LMDZ5
License: CeCILL version 2
Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
See the license file in the root directory

Also added the files defining the CeCILL version 2 license, in French
and English, at the top of the LMDZ tree.

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 3.6 KB
Line 
1module spherical
2
3  implicit none
4
5contains
6
7
8  subroutine rectsph(rect,col,azm,r,r2)
9
10    ! Converts rectangular coordinates to spherical coordinates.
11    ! (Ox) is taken as the polar axis.
12    ! Azimuth is the angle in the (Oyz) plane, measured from vector y.
13    ! If neither "r" nor "r2" is present then they are assumed to be
14    ! equal to 1.
15
16    real, intent(in):: rect(3)
17    ! (Rectangular coordinates x,y,z. Should be different from (0,0
18    ! ,0))
19
20    real, intent(out):: col ! Colatitude (in radians)
21    real, intent(out):: azm ! Azimuth (in radians). -pi < azm <= pi
22    real, intent(out), optional:: r ! Radius
23    real, intent(out), optional:: r2 ! Square of radius
24
25    ! Local variables:
26    logical p, p2
27    real d ! Radius
28    real d2 ! Square of radius
29
30    !----------------------------------------------------------------
31
32    p = present(r)
33    p2 = present(r2)
34
35    if (.not. p .and. .not. p2) then
36       ! Distance to the origin is assumed to be 1
37       col = acos(rect(1))
38    else
39       ! Compute the distance to the origin:
40       d2 = dot_product(rect,rect)
41       d = sqrt(d2)
42       col = acos(rect(1)/d)
43       if (p) r = d
44       if (p2) r2 = d2
45    end if
46
47    if (rect(2) == 0. .and. rect(3) == 0.) then
48       azm = 0.
49       ! (arbitrary value, azimuth is not well defined since vector
50       !  position is parallel
51       !  to polar axis)
52    else
53       azm = atan2(rect(3),rect(2))
54    end if
55
56  end subroutine rectsph
57
58  !*******************************************************************
59
60  real function sphbase(col,azm)
61
62    ! This function returns the matrix of the spherical base: (radial
63    ! vector, colatitude vector, azimuthal vector) in the cartesian
64    ! vector base: (x, y, z).  Vector x is assumed to be the polar
65    ! vector. Colatitude "col" and azimuth "azm" are the angular
66    ! spherical coordinates of the radial vector (azimuth is measured
67    ! from vector "y"). Note that if col = 0 (radial vector parallel
68    ! to x) then the choice of either the colatitude vector or the
69    ! azimuthal vector is arbitrary. The choice made depends on the
70    ! input value of azm. If azm is 0 then sphbase will return vector
71    ! y as the colatitude vector ("sphbase" is then the identity
72    ! matrix). Uses "sphrect", which converts spherical coordinates to
73    ! rectangular coordinates.
74
75    dimension sphbase(3,3) ! Transformation matrix (no dimension)
76
77    real, intent(in):: col ! Colatitude (in radians)
78    real, intent(in):: azm ! Azimuth (in radians)
79
80    ! Local variable:
81    real pi
82
83    !-----------------------------------------------------------------
84
85    pi = acos(-1.)
86
87    ! Radial vector:
88    sphbase(:,1) = sphrect(col, azm)
89
90    ! Colatitude vector:
91    sphbase(:,2) = sphrect(col + pi/2, azm)
92
93    ! Azimuthal vector:
94    sphbase(:,3) = (/0., - sin(azm), cos(azm)/)
95
96  end function sphbase
97
98  !*****************************************************************
99
100  real function sphrect(col,azm,r)
101
102    ! Converts spherical coordinates to rectangular coordinates. (Ox)
103    ! is taken as the polar axis. Azimuth "azm" is the angle in the
104    ! (Oyz) plane, measured from vector y.
105    ! If r is not present then we assume r = 1.
106
107    dimension sphrect(3) ! Rectangular coordinates = (x,y,z)
108
109    real, intent(in):: col ! Colatitude (in radians)
110    real, intent(in):: azm ! Azimuth (in radians)
111    real, intent(in), optional:: r ! Radius
112
113    !---------------------------------------------------------------
114
115    sphrect = (/cos(col), sin(col) * cos(azm), sin(col) * sin(azm)/)
116    if (present(r)) sphrect = sphrect * r
117
118  end function sphrect
119
120end module spherical
Note: See TracBrowser for help on using the repository browser.