⭐ 欢迎来到虫虫下载站! | 📦 资源下载 📁 资源专辑 ℹ️ 关于我们
⭐ 虫虫下载站

📄 projctns.cpp

📁 一个英国人写的GIS查看/编辑工具。支持标准的shapefile地图文件格式和coverage地图文件格式。同时可以编辑相应的dbf文件。
💻 CPP
📖 第 1 页 / 共 2 页
字号:

   // Calculate longitude
 
    X = ((1 / cos(phi)) / v);
    XI = ((1 / cos(phi)) / (6.0 * (pow3(v)))) * ((v / rho) + 2.0 * (pow2(tan(phi))));
    XII = ((1 / cos(phi)) / (120.0 * (pow5(v)))) * (5.0 + (28.0 * (pow2(tan(phi)))) + (24.0 * (pow4(tan(phi)))));
    XIIA = ((1.0 / cos(phi)) / (5040.0 * (pow7(v)))) * (61.0 + (662.0 * (pow2(tan(phi)))) + (1320.0 * (pow4(tan(phi)))) + (720 * (pow6(tan(phi)))));    
    *longitude = m_dLongitudeOrigin + ((Et * X) - ((pow3(Et)) * XI) + ((pow5(Et)) * XII) - ((pow7(Et)) * XIIA)) * 180.0 / PI;
}

////////////////////////////////////////////////////////////////////////////////
//
// Converts a string of the format DD::MM::SS NSEW to a decimal number.  The 
// flag indicates if the value is a longitude or latitude value
//

BOOL CProjection::StringAsLatLong(LPCSTR sText, double* pdValue, int iFlag)
{                                                                    
   int nDeg, nMin;
   double dDeg = 0;
   double dMin = 0, dSec = 0;
   long lSign = 1;   
   char cDir = '\0';
   char* sLatLon;
   CString s;
   
   // Retrieve values from string, nb. because the last number may be
   // a decimal any following 'e' will be interpreted as an exponent

   sLatLon = (LPSTR)sText;
   while (*sLatLon == ' ') sLatLon++;

   if (!BDNext(sLatLon, nDeg) || !BDNext(sLatLon, nMin) ||
       !BDNext(sLatLon, dSec))
   {
      sLatLon = (LPSTR)sText;
      while (*sLatLon == ' ') sLatLon++;

      if (!BDNext(sLatLon, nDeg) || !BDNext(sLatLon, dMin))
      {
        // Add support for decimal degrees on import only

         sLatLon = (LPSTR)sText;
         while (*sLatLon == ' ') sLatLon++;

         if (((CString)sLatLon).Find('.') > 0 && BDNext(sLatLon, dDeg))
         {
            nDeg = (int)dDeg;
            dMin = (dDeg - nDeg)*60;
         } else
         {
            return FALSE;
         }
      }

      nMin = (int)dMin;
      dSec = (dMin-nMin)*60;
   }         
   
   
   // Check direction

   if (BDNext(sLatLon, cDir))
   {
      if (nDeg < 0) return FALSE; // Can't have a sign and a direction

      s = cDir;
      s.MakeUpper();   
      if (s == BDString(IDS_N) && iFlag & latitude) lSign = 1;
      else if (s == BDString(IDS_S) && iFlag & latitude) lSign = -1;
      else if (s == BDString(IDS_E) && iFlag & longitude) lSign = 1;
      else if (s == BDString(IDS_W) && iFlag & longitude) lSign = -1;
      else return FALSE;
   };   

   // Validate

   if (nMin > 60.0 || dSec > 60.0) return FALSE;
   
   // Calculate decimal

   *pdValue = nDeg + nMin/60.0 + dSec/3600.0;      
   *pdValue = *pdValue * lSign;   

   // Validate range

   if (iFlag & latitude == latitude && fabs(*pdValue > 90.0)) return FALSE;
   if (iFlag & longitude == longitude && fabs(*pdValue > 180.0)) return FALSE;
   
   return TRUE;
};

///////////////////////////////////////////////////////////////////////////////

void CProjection::LatLongAsString(double deg, CString& s, int iFlag)
{
   int hour,min;
   double dSec;

   char dir;
   hour = abs((int)deg);
   min = abs((int)((fabs(deg)-hour)*60));
   dSec = fabs((((fabs(deg)-hour)*60)-min)*60);

   if (fabs(dSec - 60) < 0.000001) 
   {
      min += 1;
      dSec = 0;
   }

   if (deg >= 0 && (iFlag & latitude) == latitude) dir = 'N';
   else if (deg < 0 && (iFlag & latitude) == latitude) dir = 'S';
   else if (deg >= 0 && (iFlag & longitude) == longitude) dir = 'E';
   else if (deg < 0 && (iFlag & longitude) == longitude) dir = 'W';
   else ASSERT(FALSE);

   if (iFlag & legend)
   {
      if (iFlag & seconds) s.Format("%d'%d'%.0lf",hour,min,dSec);         
      else s.Format("%d'%d",hour,min);         
   } else
   {
      s.Format("%d'%02d\'%.1lf %c",hour,min,dSec,dir);   
   };
}

///////////////////////////////////////////////////////////////////////////////
//
// Converts a string of the format X,Y where X and Y are either latitude and
// longitude or X and Y coordinates.  If the former then a tranformation 
// is performed
//

int CProjection::StringAsCoord(CString s, CLongCoord* coordL)
{
   CCoord coord;
   if (StringAsCoord(s, &coord))
   {
      coordL->x = (long)(coord.x + 0.5); 
      coordL->y = (long)(coord.y + 0.5);   
      return TRUE;
   } else
   {
      return FALSE;
   }
}

///////////////////////////////////////////////////////////////////////////////

int CProjection::StringAsCoord(CString s, CCoord* pCoord)
{
   BOOL bOK = TRUE;
   int nRet = 0;

   CString s1, s2;   
   double dLat, dLong;
   char ch;

   // Search for comma or space

   s.TrimLeft();
   int i = s.Find(",");
   if (i == -1) i = s.Find(" ");

   if (i != -1)
   {
      s1 = s.Left(i);
      s2 = s.Mid(i+1);
      s1.TrimRight();
      s2.TrimRight();

      // Check for latitude, longitude

      if (StringAsLatLong(s1, &dLat, latitude) && 
          StringAsLatLong(s2, &dLong, longitude))
      {
         LatLonToTransMercator(dLat, dLong, &pCoord->x, &pCoord->y);         
         nRet = PR_LATLON;
      }
      else
      {
         bOK = sscanf(s1,"%lf%c",&pCoord->x, &ch) == 1 &&
               sscanf(s2,"%lf%c",&pCoord->y, &ch) == 1;
         nRet = PR_COORD;
      }
   }
   else
   {
      bOK = sscanf(s,"%lf %lf",&pCoord->x, &pCoord->y) == 2;               
   }

   if (bOK) return nRet;
   else return FALSE;
}

///////////////////////////////////////////////////////////////////////////////
//
// Converts coordinates to latitude longitude and outputs as a string
//

void CProjection::CoordAsString(CCoord coord, CString& s)
{
   double dLat, dLong;
   CString sLat, sLong;
  
   TransMercatorToLatLon(coord.x, coord.y, &dLat, &dLong);   
   LatLongAsString(dLat, sLat, latitude);
   LatLongAsString(dLong, sLong, longitude);

   s = sLat + ", " + sLong;
}

///////////////////////////////////////////////////////////////////////////////
//
// Converts a bearing in the form:
// 12N 37E 123.7 
// To a relative coordinate in metres
//

BOOL CProjection::StringAsBearing(CString sT, CLongCoord &coordL)
{
   CCoord coord;
   if (StringAsBearing(sT, coord))
   {
      coordL.x = (long)(coord.x + 0.5); 
      coordL.y = (long)(coord.y + 0.5);

      return TRUE;
   } else
   {
      return FALSE;
   }
}

///////////////////////////////////////////////////////////////////////////////

BOOL CProjection::StringAsBearing(CString s, CCoord &coord)
{
   BOOL bOK = TRUE;

   s.TrimLeft();
   int nDeg, nMin;
   double dD, dAngle;
   char cN, cE, c;
   BOOL bN = FALSE;
   BOOL bE = FALSE;

   if (sscanf(s, "%d%c%d %c%c %lf", &nDeg, &c, &nMin, &cN, &cE, &dD) == 6 && 
       (cN == 'N' || cN == 'n' || cN == 's' || cN == 'S'))
   {            
      // Determine if +ve or -ve

      if (cN == 'N' || cN == 'n')
      {
         bN = TRUE;
      }
      else if (cN == 'S' || cN == 's')
      {
         bN = FALSE;
      }
      else 
      {
         bOK = FALSE;
      }

      // Determine if +ve or -ve

      if (cE == 'E' || cE == 'e')
      {
         bE = TRUE;
      }
      else if (cE == 'W' || cE == 'w')
      {
         bE = FALSE;
      }
      else 
      {
         bOK = FALSE;
      }

      // Determine angle

      dAngle = nDeg + nMin/60;
      if (!bN && bE) dAngle = 90 + 90 - dAngle;
      if (!bN && !bE) dAngle = 180 + dAngle;
      if (bN && !bE) dAngle = 270 + 90 - dAngle;
      dAngle = dAngle/90.0 * PI/2;

      // Determine bearing

      coord.y = cos(dAngle) * dD;
      coord.x = sin(dAngle) * dD;
           
   } else
   {
      bOK = FALSE;
   }

   return bOK;
}

///////////////////////////////////////////////////////////////////////////////

CString CProjection::GetProjectionName()
{
   return this->m_sName;
}

CBDProjection& CProjection::GetProjection()
{
   return (CBDProjection&)*this;
}

///////////////////////////////////////////////////////////////////////////////
//
// Return distance between coordinates in metres
//
// www.mathforum.org
//
// cos(AOB) = cos(dLatA)cos(dLatB)cos(lonB-lonA)+sin(dLatA)sin(dLatB)
//
// This gives AOB, and the great circle distance between A and
// B will be
//               R*(AOB)     with AOB in radians.
//
// NOTE: Treats earth as sphere not elipsoid!
//

int CProjection::GetDistance(CLongCoord coord1, CLongCoord coord2)
{
   // For tranverse mercator this is simple as in metres

  if (!IsLatLon())
  {     
     // Force use of floating point arithmetic to remove rounding errors

     return (int)sqrt(square(coord2.x - coord1.x) + square(coord2.y - coord1.y));
  } 
  
  // Calculate for lat/long coordinates

  else
  {  
    // Convert to latitude/longitude
     
     double dLat1, dLon1, dLat2, dLon2;
     TransMercatorToLatLon(coord1.x, coord1.y, &dLat1, &dLon1);
     TransMercatorToLatLon(coord2.x, coord2.y, &dLat2, &dLon2);

     // Convert to radians

     dLat1 = dLat1 / 180 * PI;
     dLon1 = dLon1 / 180 * PI;
     dLat2 = dLat2 / 180 * PI;
     dLon2 = dLon2 / 180 * PI;

     double c = acos( (cos(dLon1-dLon2)-1)*cos(dLat1)*cos(dLat2) + cos(dLat1-dLat2) );

     // Determine radius for this latitude

     double r = EARTHRADIUS * (1 - 0.0033493 * square(sin((dLat1+dLat2)/2)));
          
     return (int)(r * c);

  }
}

⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -