Quantum GIS API Documentation  1.7.4
src/core/qgsscalecalculator.cpp
Go to the documentation of this file.
00001 /***************************************************************************
00002                               qgsscalecalculator.h
00003                  Calculates scale based on map extent and units
00004                               -------------------
00005   begin                : May 18, 2004
00006   copyright            : (C) 2004 by Gary E.Sherman
00007   email                : sherman at mrcc.com
00008  ***************************************************************************/
00009 
00010 /***************************************************************************
00011  *                                                                         *
00012  *   This program is free software; you can redistribute it and/or modify  *
00013  *   it under the terms of the GNU General Public License as published by  *
00014  *   the Free Software Foundation; either version 2 of the License, or     *
00015  *   (at your option) any later version.                                   *
00016  *                                                                         *
00017  ***************************************************************************/
00018 /* $Id$ */
00019 
00020 #include <cmath>
00021 #include "qgslogger.h"
00022 #include "qgsrectangle.h"
00023 #include "qgsscalecalculator.h"
00024 
00025 QgsScaleCalculator::QgsScaleCalculator( double dpi, QGis::UnitType mapUnits )
00026     : mDpi( dpi ), mMapUnits( mapUnits )
00027 {}
00028 
00029 QgsScaleCalculator::~QgsScaleCalculator()
00030 {}
00031 
00032 void QgsScaleCalculator::setDpi( double dpi )
00033 {
00034   mDpi = dpi;
00035 }
00036 double QgsScaleCalculator::dpi()
00037 {
00038   return mDpi;
00039 }
00040 
00041 void QgsScaleCalculator::setMapUnits( QGis::UnitType mapUnits )
00042 {
00043   QgsDebugMsg( QString( "Map units set to %1" ).arg( QString::number( mapUnits ) ) );
00044   mMapUnits = mapUnits;
00045 }
00046 
00047 QGis::UnitType QgsScaleCalculator::mapUnits() const
00048 {
00049   QgsDebugMsgLevel( QString( "Map units returned as %1" ).arg( QString::number( mMapUnits ) ), 3 );
00050   return mMapUnits;
00051 }
00052 
00053 double QgsScaleCalculator::calculate( const QgsRectangle &mapExtent, int canvasWidth )
00054 {
00055   double conversionFactor = 0;
00056   double delta = 0;
00057   // calculation is based on the map units and extent, the dpi of the
00058   // users display, and the canvas width
00059   switch ( mMapUnits )
00060   {
00061     case QGis::Meters:
00062       // convert meters to inches
00063       conversionFactor = 39.3700787;
00064       delta = mapExtent.xMaximum() - mapExtent.xMinimum();
00065       break;
00066     case QGis::Feet:
00067       conversionFactor = 12.0;
00068       delta = mapExtent.xMaximum() - mapExtent.xMinimum();
00069       break;
00070     case QGis::DecimalDegrees:
00071       // degrees require conversion to meters first
00072       conversionFactor = 39.3700787;
00073       delta = calculateGeographicDistance( mapExtent );
00074       break;
00075     case QGis::DegreesMinutesSeconds:
00076       // degrees require conversion to meters first
00077       conversionFactor = 39.3700787;
00078       delta = calculateGeographicDistance( mapExtent );
00079       break;
00080     case QGis::DegreesDecimalMinutes:
00081       // degrees require conversion to meters first
00082       conversionFactor = 39.3700787;
00083       delta = calculateGeographicDistance( mapExtent );
00084       break;
00085     default:
00086       Q_ASSERT( "bad map units" );
00087       break;
00088   }
00089   QgsDebugMsg( "Using conversionFactor of " + QString::number( conversionFactor ) );
00090   if ( canvasWidth == 0 || mDpi == 0 )
00091   {
00092     QgsDebugMsg( "Can't calculate scale from the input values" );
00093     return 0;
00094   }
00095   double scale = ( delta * conversionFactor ) / (( double )canvasWidth / mDpi );
00096   return scale;
00097 }
00098 
00099 
00100 double  QgsScaleCalculator::calculateGeographicDistance( const QgsRectangle &mapExtent )
00101 {
00102   // need to calculate the x distance in meters
00103   // We'll use the middle latitude for the calculation
00104   // Note this is an approximation (although very close) but calculating scale
00105   // for geographic data over large extents is quasi-meaningless
00106 
00107   // The distance between two points on a sphere can be estimated
00108   // using the Haversine formula. This gives the shortest distance
00109   // between two points on the sphere. However, what we're after is
00110   // the distance from the left of the given extent and the right of
00111   // it. This is not necessarily the shortest distance between two
00112   // points on a sphere.
00113   //
00114   // The code below uses the Haversine formula, but with some changes
00115   // to cope with the above problem, and also to deal with the extent
00116   // possibly extending beyond +/-180 degrees:
00117   //
00118   // - Use the Halversine formula to calculate the distance from -90 to
00119   //   +90 degrees at the mean latitude.
00120   // - Scale this distance by the number of degrees between
00121   //   mapExtent.xMinimum() and mapExtent.xMaximum();
00122   // - For a slight improvemnt, allow for the ellipsoid shape of earth.
00123 
00124 
00125   // For a longitude change of 180 degrees
00126   double lat = ( mapExtent.yMaximum() + mapExtent.yMinimum() ) * 0.5;
00127   const static double rads = ( 4.0 * atan( 1.0 ) ) / 180.0;
00128   double a = pow( cos( lat * rads ), 2 );
00129   double c = 2.0 * atan2( sqrt( a ), sqrt( 1.0 - a ) );
00130   const static double ra = 6378000; // [m]
00131   // The eccentricity. This comes from sqrt(1.0 - rb*rb/(ra*ra)) with rb set
00132   // to 6357000 m.
00133   const static double e = 0.0810820288;
00134   double radius = ra * ( 1.0 - e * e ) /
00135                   pow( 1.0 - e * e * sin( lat * rads ) * sin( lat * rads ), 1.5 );
00136   double meters = ( mapExtent.xMaximum() - mapExtent.xMinimum() ) / 180.0 * radius * c;
00137 
00138   QgsDebugMsg( "Distance across map extent (m): " + QString::number( meters ) );
00139 
00140   return meters;
00141 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines