Convertir les coordonnées OSGB36 en Latitude / Longitude

Question

Je souhaite convertir les coordonnées OSGB 36 en WGS 84 (c.-à-d. &" standard & "; latitude et longitude), afin de les tracer dans un fichier KML pour Google Earth.

Quelle serait la meilleure façon de s'y prendre? J'implémente en VB .NET.

Je devrais probablement ajouter que ma question n'est pas & "Comment puis-je écrire un fichier KML? &"; Ma question est & «Comment convertir les deux systèmes de coordonnées? &« !!

J'espérais qu'il y aurait une bibliothèque que je pourrais utiliser plutôt que de créer ma propre fonction - il me semble que ce serait le genre de chose que quelqu'un d'autre aurait mise en oeuvre.

Answer 1

La société pour laquelle je travaille vient d'ouvrir une bibliothèque à partir de cette source: http: // code.google.com/p/geocoordconversion/

Answer 2

Vous devez implémenter une transformation Helmert . J'ai écrit une conversion en Javascript que vous pourrez peut-être adapter. L'algorithme utilisé par le script pour les conversions WGS84-OSGB36 est dérivé d'un tableur OSGB avec autorisation. La précision de conversion est de l’ordre de 7 m pour 90% de la Grande-Bretagne et devrait être similaire à la conversion effectuée par un récepteur GPS typique.

Voir la documentation et source pour plus de détails.

Edit: vous voudrez peut-être consulter ce OCX , qui comprend source.

Answer 3

Premièrement, selon this Lien vers cette page à partir de OSGB 36 :

  

Mythe 4: & # 8216; Il existe des formules mathématiques exactes permettant de changer de système de coordonnées & # 8217;

Deuxièmement, à partir du même lien: & "; D'un système de coordonnées à un autre: transformations géodésiques " inclut une section & "; Conversion approximative de WGS84 à OSGB36 / ODN &";

EDIT: Notez que lorsque le système d'exploitation dit & "approximativement &", ils entendent par des erreurs > 5 m.

Answer 4

//=======================================================================
/* Project: Geocode.Service
*  Author: Steve Loughran
*  Copyright (C) 2000 Steve Loughran
*  See license.txt or license.html for license and copyright information 
*  RCS $Header: /cvsroot/geocode/geocode.net/src/library/Osgb.cs,v 1.4 2002/04/23 05:12:37 steve_l Exp $
*  jedit:mode=csharp:tabSize=4:indentSize=4:syntax=true:
*/
//=======================================================================


using System;

namespace net.sourceforge.geocode {

/// <summary>
/// OSGB conversion routines. xlated from C++ to Java to C#
/// </summary>

public class OsgbGridReference: GeoMath
{

private string _gridSquare="";
private long _easting=0;
private long _northing=0;

/// <summary>
/// empty constructor
/// </summary>
public OsgbGridReference() {}

/// <summary>
/// constructor from gridref
/// </summary>
public OsgbGridReference(String gridSquare, 
  long easting,
  long northing) {
 SetGridRef(gridSquare,northing,easting);
}

/// <summary>
/// constructor from co-ordinates
/// </summary>
public OsgbGridReference(double latitude, double longitude) {
 SetPosition(latitude,longitude);
}

/// <summary>
/// constructor from position
/// </summary>
public OsgbGridReference(Position position)
 : this(position.Latitude,position.Longitude) {
}

/// <summary>grid square property</summary>    
public string GridSquare {
 get {return _gridSquare;}
 set {_gridSquare=value;}
}
/// <summary>northing property</summary>    

public long Northing {
 get {return _northing;}
 set {_northing=value;} 
}

/// <summary>easting property</summary>    
public long Easting {
 get {return _easting;}
 set {_easting=value;} 
}

/// <summary>
/// set the grid refernce
/// </summary>
/// <returns> </returns>

public void SetGridRef(String gridSquare, 
  long northing, 
  long easting) {
 _gridSquare=gridSquare;
 _northing=northing;
 _easting=easting;
}

/// <summary>
///  rudimentary validity test. A better one is to
/// extract the position
/// </summary>
/// <returns>true iff there is a gridsquare </returns>

public bool IsValid() 
 {return _gridSquare.Length==2;}


/// <summary>
/// get a position object from a location
/// </summary>
/// <returns> Position </returns>

public Position ToPosition() {
 double lat,lon;
 ConvertOSGBtoLL(_gridSquare,_northing,_easting,
   out lat, out lon);
 Position p=new Position(lat,lon);
 p.Name=ToString();
 return p;
}   

/// <summary>
/// set a gridref from a lat/long tuple
/// </summary>
/// <returns> success flag </returns>

public bool SetPosition(double latitude, double longitude) {
 _gridSquare=ConvertLLtoOSGB(latitude,
   longitude,
   out _northing,
   out _easting);
 return IsValid();
}

/// <summary>
/// set a gridref from a position
/// </summary>
/// <returns> success flag </returns>

public bool SetPosition(Position position) {
 return SetPosition(position.Latitude,position.Longitude);
}

/// <summary>
///  stringify
/// </summary>

public override string ToString() {
 return String.Format("{0} {1:000} {2:000}",
  _gridSquare,Easting,Northing);
}

/// <summary>
///  equality test: works on lat and long
/// </summary>

public override bool Equals(object o) {
 OsgbGridReference pos=(OsgbGridReference)o;
 return _gridSquare==pos._gridSquare &&
  _northing==pos._northing &&
  _easting==pos._easting;
}

/// <summary>
/// hash code builder
/// </summary>

public override int GetHashCode() {
        return (int)(_easting+_northing); 
}

/// <summary>
/// determines base co-ordinates of a square like "ST"
/// </summary>
///    <parameter name="OSGBGridSquare"> square </parameter>
///    <parameter name="easting"> easting</parameter>  
///    <parameter name="northing"> northing</parameter> 
/// <returns>true if the coordinates were in range</returns>

static bool ConvertOSGBSquareToRefCoords(string OSGBGridSquare,
         out long  easting,
                           out long  northing) {
 int pos, x_multiplier=0, y_multiplier=0;
 string GridSquare = "VWXYZQRSTULMNOPFGHJKABCDE";
 bool trouble=false;
 long east,north;
 easting=northing=0;
 //find 500km offset
    OSGBGridSquare=OSGBGridSquare.ToUpper();
 char ch = OSGBGridSquare[0];
 switch(ch) {
  case 'S': x_multiplier = 0; y_multiplier = 0; break;
  case 'T': x_multiplier = 1; y_multiplier = 0; break;
  case 'N': x_multiplier = 0; y_multiplier = 1; break;
  case 'O': x_multiplier = 1; y_multiplier = 1; break;
  case 'H': x_multiplier = 0; y_multiplier = 2; break;
  case 'J': x_multiplier = 1; y_multiplier = 2; break;
        default: trouble=true; break;
 }
    if(!trouble) {
  east=x_multiplier * 500000L;
  north=y_multiplier * 500000L;
  //find 100km offset and add to 500km offset to get coordinate of
  //square point is in
  char subsquare=OSGBGridSquare[1];
  pos = GridSquare.IndexOf(subsquare);
     if(pos>-1) {
   east += ((pos % 5) * 100000L);
   north += ((pos / 5) * 100000L);
     easting=east;
   northing=north;
     }
     else {
         trouble=true;
     }
    }
    return !trouble;
}

///<summary>
///convert a internal OSGB coord to lat/long
///Equations from USGS Bulletin 1532
///East Longitudes are positive, West longitudes are negative.
///North latitudes are positive, South latitudes are negative
///Lat and Long are in decimal degrees.
///Written by Chuck Gantz- chuck.gantz@globalstar.com
///</summary>
///    <parameter name="OSGBEasting">easting </parameter> 
///   <parameter name="OSGBEasting">northing </parameter> 
///   <parameter name="OSGBZone"> gridsquare</parameter>  
///   <parameter name="latitude"> latitude</parameter> 
///   <parameter name="longitude"> longitude</parameter> 

static void ConvertOSGBtoLL(string OSGBZone,
   double OSGBNorthing,
   double OSGBEasting,
   out double latitude,
   out double longitude) {
 double k0 = 0.9996012717;
 double a;
 double eccPrimeSquared;
 double N1, T1, C1, R1, D, M;
 double LongOrigin = -2;
 double LatOrigin = 49;
 double LatOriginRad = LatOrigin * deg2rad;
 double mu, phi1, phi1Rad;
 double x, y;
    long northing;
    long easting;

 //Airy model of the ellipsoid.
 double majoraxis = a = 6377563.396;
 double minoraxis = 6356256.91;
 double eccSquared = (majoraxis * majoraxis - minoraxis * minoraxis) /
       (majoraxis * majoraxis);
 double e1 = (1-sqrt(1-eccSquared))/(1+sqrt(1-eccSquared));
 //only calculate M0 once since it is based on the origin of the OSGB projection, which is fixed
 double  M0 = a*((1 - eccSquared/4 - 3*eccSquared*eccSquared/64 - 5*eccSquared*eccSquared*eccSquared/256)*LatOriginRad
    - (3*eccSquared/8 + 3*eccSquared*eccSquared/32 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(2*LatOriginRad)
    + (15*eccSquared*eccSquared/256 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(4*LatOriginRad)
    - (35*eccSquared*eccSquared*eccSquared/3072)*sin(6*LatOriginRad));
 ConvertOSGBSquareToRefCoords(OSGBZone, out easting, out northing);
 //remove 400,000 meter false easting for longitude
 x = OSGBEasting - 400000.0 + easting;
 //remove 100,000 meter false easting for longitude
 y = OSGBNorthing + 100000.0 + northing;
 eccPrimeSquared = (eccSquared)/(1-eccSquared);
 M = M0 + y / k0;
 mu = M/(a*(1-eccSquared/4-3*eccSquared*eccSquared/
      64-5*eccSquared*eccSquared*eccSquared/256));
 phi1Rad = mu + (3*e1/2-27*e1*e1*e1/32)*sin(2*mu)
    + (21*e1*e1/16-55*e1*e1*e1*e1/32)*sin(4*mu)
    +(151*e1*e1*e1/96)*sin(6*mu);
 phi1 = phi1Rad*rad2deg;
 N1 = a/sqrt(1-eccSquared*sin(phi1Rad)*sin(phi1Rad));
 T1 = tan(phi1Rad)*tan(phi1Rad);
 C1 = eccPrimeSquared*cos(phi1Rad)*cos(phi1Rad);
 R1 = a*(1-eccSquared)/pow(1-eccSquared*sin(phi1Rad)*sin(phi1Rad), 1.5);
 D = x/(N1*k0);
 latitude = phi1Rad - (N1*tan(phi1Rad)/R1)*(D*D/2-(5+3*T1+10*C1-4*C1*C1-9*eccPrimeSquared)*D*D*D*D/24
     +(61+90*T1+298*C1+45*T1*T1-252*eccPrimeSquared-3*C1*C1)*D*D*D*D*D*D/720);
 latitude *= rad2deg;
 longitude = (D-(1+2*T1+C1)*D*D*D/6+(5-2*C1+28*T1-3*C1*C1+8*eccPrimeSquared+24*T1*T1)
     *D*D*D*D*D/120)/cos(phi1Rad);
 longitude = LongOrigin + longitude * rad2deg;
}


/// <summary>
/// converts lat/long to OSGB coords.  Equations from USGS Bulletin 1532
/// East Longitudes are positive, West longitudes are negative.
/// North latitudes are positive, South latitudes are negative
/// Lat and Long are in decimal degrees
/// </summary>
///  Written by Chuck Gantz- chuck.gantz@globalstar.com
///   <parameter name="latitude"> IN latitude</parameter> 
///   <parameter name="longitude">IN longitude </parameter> 
///   <parameter name="OSGBEasting"> OUT easting</parameter> 
///  <parameter name="OSGBNorthing"> OUT northing</parameter> 

static public string ConvertLLtoOSGB(double latitude,
    double longitude,                
                out long OSGBNorthing,
    out long OSGBEasting) {
 double a;
 double eccSquared;
 double k0 = 0.9996012717;
 double LongOrigin = -2;
 double LongOriginRad = LongOrigin * deg2rad;
 double LatOrigin = 49;
 double LatOriginRad = LatOrigin * deg2rad;
 double eccPrimeSquared;
 double N, T, C, A, M;
 double LatRad = latitude*deg2rad;
 double LongRad = longitude*deg2rad;
 double easting, northing;
 double majoraxis = a = 6377563.396;//Airy
 double minoraxis = 6356256.91;//Airy
 eccSquared = (majoraxis * majoraxis - minoraxis * minoraxis) /
       (majoraxis * majoraxis);
 //only calculate M0 once since it is based on the origin
 //of the OSGB projection, which is fixed
 double  M0 = a*((1  - eccSquared/4 - 3*eccSquared*eccSquared/64 - 5*eccSquared*eccSquared*eccSquared/256)*LatOriginRad
    - (3*eccSquared/8 + 3*eccSquared*eccSquared/32 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(2*LatOriginRad)
    + (15*eccSquared*eccSquared/256 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(4*LatOriginRad)
    - (35*eccSquared*eccSquared*eccSquared/3072)*sin(6*LatOriginRad));
 eccPrimeSquared = (eccSquared)/(1-eccSquared);
 N = a/sqrt(1-eccSquared*sin(LatRad)*sin(LatRad));
 T = tan(LatRad)*tan(LatRad);
 C = eccPrimeSquared*cos(LatRad)*cos(LatRad);
 A = cos(LatRad)*(LongRad-LongOriginRad);
 M = a*((1 - eccSquared/4 - 3*eccSquared*eccSquared/64- 5*eccSquared*eccSquared*eccSquared/256)*LatRad
    - (3*eccSquared/8 + 3*eccSquared*eccSquared/32 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(2*LatRad)
    + (15*eccSquared*eccSquared/256 + 45*eccSquared*eccSquared*eccSquared/1024)*sin(4*LatRad)
    - (35*eccSquared*eccSquared*eccSquared/3072)*sin(6*LatRad));
 easting = (double)(k0*N*(A+(1-T+C)*A*A*A/6
    + (5-18*T+T*T+72*C-58*eccPrimeSquared)*A*A*A*A*A/120));
 easting += 400000.0; //false easting
 northing = (double)(k0*(M-M0+N*tan(LatRad)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24
     + (61-58*T+T*T+600*C-330*eccPrimeSquared)*A*A*A*A*A*A/720)));
 northing -= 100000.0;//false northing
 return ConvertCoordsToOSGBSquare(easting, northing,out OSGBEasting, out OSGBNorthing);
}


/// <summary>
/// convert a internal OSGB coord to a gridsquare and internal values.
/// </summary>
///    <parameter name="easting"> easting</parameter> 
///    <parameter name="northing"> northing</parameter>
///    <parameter name="OSGBEasting">OSGBEasting out </parameter>
///    <parameter name="OSGBNorthing">OSGBNorthing out </parameter>

static string ConvertCoordsToOSGBSquare(double easting,
     double northing,
                    out long  OSGBEasting,
                    out long  OSGBNorthing)
{
 string GridSquare = "VWXYZQRSTULMNOPFGHJKABCDE";
 long posx, posy; //positions in grid
 OSGBEasting = (long)(easting + 0.5); //round to nearest int
 OSGBNorthing = (long)(northing + 0.5); //round to nearest int
 string OSGBGridSquare="";

 //find correct 500km square
 posx = OSGBEasting / 500000L;
 posy = OSGBNorthing / 500000L;
 if(posx<0 || posx>4 || posy<0 || posy>4)
  return "";
 long offset=posx + posy * 5 + 7;
    OSGBGridSquare+= GridSquare[(int)offset];

 //find correct 100km square
 posx = OSGBEasting % 500000L;//remove 500 km square
 posy = OSGBNorthing % 500000L;//remove 500 km square
 posx = posx / 100000L;//find 100 km square
 posy = posy / 100000L;//find 100 km square
 if(posx<0 || posx>4 || posy<0 || posy>4)
  return "";
 offset=posx + posy * 5; 
 OSGBGridSquare+= GridSquare[(int)offset];

 //remainder is northing and easting
 OSGBNorthing = OSGBNorthing % 500000L;
 OSGBNorthing = OSGBNorthing % 100000L;
 OSGBEasting = OSGBEasting % 500000L;
 OSGBEasting = OSGBEasting % 100000L;
    return OSGBGridSquare;
}



/// <summary>
/// turn the latitude and longitude into a string
/// </summary>
///    <parameter name="latitude"> lat</parameter>
///    <parameter name="longitude"> long</parameter> 
///    <parameter name="infill"> text between coordinates</parameter>  
///    <returns>return something like E004 N123</returns>

static string GetSensibleLatLongstring(double latitude,
  double longitude,
        int decimals,
        string infill) {
    bool bNorth=latitude>=0;
 bool bWest=longitude<=0;
    latitude=Math.Abs(latitude);
    longitude=Math.Abs(longitude);
    double multiplier=(int)pow(10,decimals);
 int hiLat=(int)latitude;
    int lowLat=(int)((latitude-hiLat)*multiplier);
    double newLat=lowLat/multiplier+hiLat;
 int hiLong=(int)longitude;
    int lowLong=(int)((longitude-hiLong)*multiplier);
    double newLong=lowLong/multiplier+hiLong;
 return (bNorth?"N":"S")+newLat+infill+
  (bWest?"W":"E")+newLong;
}




/* legacy java test harness
  public static void main(string[] args)
 {
    string message;
    if(args.length!=3)
     {
        message="need a grid reference like ST 767 870";
        }
    else
     {
        LongRef north=new LongRef();
        LongRef east=new LongRef();
        bool b=ConvertOSGBSquareToRefCoords(args[0],east,north);
        double easting=Double.valueOf(args[1]).doubleValue();
        double northing=Double.valueOf(args[2]).doubleValue();
        DoubleRef rlatitude=new DoubleRef ();
        DoubleRef rlongitude=new DoubleRef ();
        OSGBtoLL(easting,northing,args[0],rlatitude,rlongitude);
        double latitude=rlatitude.getValue();
        double longitude=rlongitude.getValue();

        bool bNorth=latitude>=0;
        bool bWest=longitude<=0;

  message="Latitude & Longitude ="+latitude+" , " + longitude;
        message+="\r\n or "+GetSensibleLatLongstring(latitude,
             longitude,
                            3,
                            " ");
  string square=new string();
  square=LLtoOSGB(latitude,longitude,north,east);
  message+="\r\n Grid ="+square+" "+east+" "+north;
//       message="evaluation failure on "+args[0];
        }
        System.out.print(message);
    }
*/


}; //class
};//namespace

Answer 5

Nous utilisons la bibliothèque proj.4 pour WGS84 < - > ; OSGB36 & Lt; - & Gt; Les transformations de coordonnées OSGBGRID et cela fonctionne très bien. Mais nous utilisons le C ++, donc je ne sais pas si vous pouvez le faire fonctionner sous VB.NET. Il peut y avoir des wrappers ou quelque chose? (Sur le lien ci-dessus, il est fait mention de PROJ.4 VB Wrappers).