7 !! Computation of position (x,y) in the numerical domain for longitude lambda
8 !! and latitude phi (in rad).
12 !! Copyright 2009-2013 Ralf Greve, Reinhard Calov, Alex Robinson
16 !! This file is part of SICOPOLIS.
18 !! SICOPOLIS is free software: you can redistribute it and/or modify
19 !! it under the terms of the GNU General Public License as published by
20 !! the Free Software Foundation, either version 3 of the License, or
21 !! (at your option) any later version.
23 !! SICOPOLIS is distributed in the hope that it will be useful,
24 !! but WITHOUT ANY WARRANTY; without even the implied warranty of
25 !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
26 !! GNU General Public License for more details.
28 !! You should have received a copy of the GNU General Public License
29 !! along with SICOPOLIS. If not, see <http://www.gnu.org/licenses/>.
31 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
34 !> Computation of position (x,y) in the numerical domain for longitude lambda
35 !! and latitude phi (in rad).
36 !<------------------------------------------------------------------------------
37 subroutine num_coord(lambda_val, phi_val, x_val, y_val)
44 real(dp),
intent(in) :: lambda_val, phi_val
46 real(dp),
intent(out) :: x_val, y_val
48 #if (GRID==0 || GRID==1)
53 lambda0, phi0, x_val, y_val)
67 !> Forward stereographic projection for an ellipsoidal planet.
68 !<------------------------------------------------------------------------------
70 lambda0, phi0, x_val, y_val)
75 real(dp),
intent(in) :: lambda_val, phi_val, a, b, lambda0, phi0
76 real(dp),
intent(out) :: x_val, y_val
79 integer(i4b) :: sign_phi0
80 real(dp) :: phi_aux, phi0_aux
81 real(dp) :: e, mc, t, tc, kp, rho, phi_p
82 real(dp) :: sinphi0, sinlambda0, cosphi0, coslambda0
84 real(dp),
parameter :: pi = 3.141592653589793_dp, &
89 else if (phi0 < (-eps))
then
92 stop
' num_coord: phi0 must be different from zero!'
95 phi_aux = phi_val * sign_phi0
96 phi0_aux = phi0 * sign_phi0
98 e=sqrt((a**2-b**2)/(a**2))
100 sinphi0 = sin(phi0_aux)
101 sinlambda0 = sin(lambda0)
102 cosphi0 = cos(phi0_aux)
103 coslambda0 = cos(lambda0)
105 mc=cosphi0/sqrt(1.0_dp-e*e*sinphi0*sinphi0)
106 t=sqrt(((1.0_dp-sin(phi_aux))/(1.0_dp+sin(phi_aux)))* &
107 ((1.0_dp+e*sin(phi_aux))/(1.0_dp-e*sin(phi_aux)))**e)
108 tc=sqrt(((1.0_dp-sinphi0)/(1.0_dp+sinphi0))* &
109 ((1.0_dp+e*sinphi0)/(1.0_dp-e*sinphi0))**e)
112 x_val = rho*sin(lambda_val-lambda0)
113 y_val = -sign_phi0 * rho*cos(lambda_val-lambda0)