! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ ! + + ! + glimmer_map_types.f90 - part of the Glimmer-CISM ice model + ! + + ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ ! ! Copyright (C) 2004, 2005, 2006, 2007, 2008, 2009, 2010 ! Glimmer-CISM contributors - see AUTHORS file for list of contributors ! ! This file is part of Glimmer-CISM. ! ! Glimmer-CISM is free software: you can redistribute it and/or modify ! it under the terms of the GNU General Public License as published by ! the Free Software Foundation, either version 2 of the License, or (at ! your option) any later version. ! ! Glimmer-CISM is distributed in the hope that it will be useful, ! but WITHOUT ANY WARRANTY; without even the implied warranty of ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the ! GNU General Public License for more details. ! ! You should have received a copy of the GNU General Public License ! along with Glimmer-CISM. If not, see . ! ! Glimmer-CISM is hosted on BerliOS.de: ! https://developer.berlios.de/projects/glimmer-cism/ ! ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ #ifdef HAVE_CONFIG_H #include "config.inc" #endif !> This module contains derived types. !! !! \todo move math constants to constants module module glimmer_map_types use glimmer_global, only: rk implicit none !> derived type holding all know map projections. This simulates inheritance type glimmap_proj logical :: found = .false. type(proj_laea), pointer :: laea => NULL() !< Pointer to Lambert azimuthal equal area type type(proj_aea), pointer :: aea => NULL() !< Pointer to Albers equal area conic type type(proj_lcc), pointer :: lcc => NULL() !< Pointer to Lambert conic conformal type type(proj_stere), pointer :: stere => NULL() !< Pointer to Stereographic type end type glimmap_proj !------------------------------------------------------------- !> Lambert Azimuthal Equal Area type proj_laea real(rk) :: longitude_of_central_meridian !< longitude of central meridian real(rk) :: latitude_of_projection_origin !< latitude of projection origin real(rk) :: false_easting !< false easting real(rk) :: false_northing !< false northing real(rk) :: sinp !< Sine of latitude_of_projection_origin real(rk) :: cosp !< Cosine of latitude_of_projection_origin integer :: pole !< Set to 1 for N pole, -1 for S pole, 0 otherwise end type proj_laea !------------------------------------------------------------- !> Albers Equal-Area Conic type proj_aea real(rk),dimension(2) :: standard_parallel !< two standard parallels real(rk) :: longitude_of_central_meridian !< longitude of central meridian real(rk) :: latitude_of_projection_origin !< latitude of projection origin real(rk) :: false_easting !< false easting real(rk) :: false_northing !< false northing real(rk) :: rho0 !< Convenience constant real(rk) :: rho0_R !< Convenience constant (is rho0/EQ_RAD) real(rk) :: c !< Convenience constant real(rk) :: n !< Convenience constant real(rk) :: i_n !< Convenience constant (inverse of n) end type proj_aea !------------------------------------------------------------- !> Lambert Conic Conformal type proj_lcc real(rk),dimension(2) :: standard_parallel !< two standard parallels real(rk) :: longitude_of_central_meridian !< longitude of central meridian real(rk) :: latitude_of_projection_origin !< latitude of projection origin real(rk) :: false_easting !< false easting real(rk) :: false_northing !< false northing real(rk) :: rho0 !< Convenience constant real(rk) :: f !< Convenience constant real(rk) :: n !< Convenience constant real(rk) :: i_n !< Convenience constant (inverse of n) end type proj_lcc !------------------------------------------------------------- !> Stereographic projection derived type type proj_stere real(rk) :: longitude_of_central_meridian !< longitude of central meridian real(rk) :: latitude_of_projection_origin !< latitude of projection origin real(rk) :: scale_factor_at_proj_origin = 0. !< scale factor at origin real(rk) :: standard_parallel = 0. !< a standard parallel real(rk) :: false_easting !< false easting real(rk) :: false_northing !< false northing integer :: pole !< Set to 1 for N pole, -1 for S pole, 0 otherwise logical :: equatorial !< Set true if equatorial aspect real(rk) :: k0 !< scale factor or std par converted to scale factor real(rk) :: ik0 !< inverse of k0 real(rk) :: sinp !< sin of latitude_of_projection_origin real(rk) :: cosp !< cos of latitude_of_projection_origin end type proj_stere ! Global mapping parameters ---------------------------------- real(rk),parameter :: pi = 3.141592654 !< The value of $\pi$. real(rk),parameter :: M_PI_4 = pi/4 !< The value of $\pi/4$. real(rk),parameter :: M_PI_2 = pi/2 !< The value of $\pi/2$. real(rk),parameter :: D2R = pi/180.0 !< Degrees-to-radians conversion factor. real(rk),parameter :: R2D = 180.0/pi !< Radians-to-degrees conversion factor. real(rk),parameter :: EQ_RAD = 6.37e6 !< Radius of the earth (m) real(rk),parameter :: i_EQ_RAD = 1.0_rk/EQ_RAD !< Inverse radius of the earth (m^-1) real(rk),parameter :: CONV_LIMIT = 1.0e-8 !< Convergence limit (a small number). integer, parameter :: GMAP_LAEA=1 !< ID for Lambert azimuthal equal area projection integer, parameter :: GMAP_AEA=2 !< ID for Lambert azimuthal equal area projection integer, parameter :: GMAP_LCC=3 !< ID for Lambert conformal conic projection integer, parameter :: GMAP_STERE=4 !< ID for stereographic projection contains !> return true if structure contains a known projection function glimmap_allocated(proj) implicit none type(glimmap_proj) :: proj logical glimmap_allocated glimmap_allocated = proj%found end function glimmap_allocated !> This is incomplete diagnostics code to output full !! content of projection type. Only does !! Stereographic projections so far. subroutine glimmap_diag(proj) use glimmer_log type(glimmap_proj) :: proj if (associated(proj%stere)) then call glimmap_diag_stere(proj%stere) else call write_log('Stereographic projection not found') end if end subroutine glimmap_diag !> print out parameters of Stereographic projection subroutine glimmap_diag_stere(params) use glimmer_log use glimmer_global, only : msg_length type(proj_stere) :: params character(len=msg_length) :: message call write_log('***** Stereographic *****') write(message,*)'longitude_of_central_meridian:', params%longitude_of_central_meridian call write_log(message) write(message,*)'latitude_of_projection_origin:', params%latitude_of_projection_origin call write_log(message) write(message,*)'scale_factor_at_proj_origin:', params%scale_factor_at_proj_origin call write_log(message) write(message,*)'standard_parallel:', params%standard_parallel call write_log(message) write(message,*)'false_easting:', params%false_easting call write_log(message) write(message,*)'false_northing:', params%false_northing call write_log(message) write(message,*)'pole:', params%pole call write_log(message) write(message,*)'equatorial:', params%equatorial call write_log(message) write(message,*)'k0:', params%k0 call write_log(message) write(message,*)'ik0:', params%ik0 call write_log(message) write(message,*)'sinp:', params%sinp call write_log(message) write(message,*)'cosp:', params%cosp call write_log(message) end subroutine glimmap_diag_stere end module glimmer_map_types