Needed to convert LAMBERT 1972 (EPSG:31370) coordinates to WGS84 (EPSG:4326) for a project I’m working on in order to use the coordinates with Google Maps *(which uses a SphericalMercator (EPSG:900913) projection)*.

Took me quite a while to find the needed algorithm, so hereby I’m listing them (and am providing the JavaScript implementation):

JavaScript:

```
var lambert72toWGS84 = function(x, y){
var newLongitude, newLatitude;
var n = 0.77164219,
F = 1.81329763,
thetaFudge = 0.00014204,
e = 0.08199189,
a = 6378388,
xDiff = 149910,
yDiff = 5400150,
theta0 = 0.07604294;
var xReal = xDiff - x,
yReal = yDiff - y;
var rho = Math.sqrt(xReal * xReal + yReal * yReal),
theta = Math.atan(xReal / -yReal);
newLongitude = (theta0 + (theta + thetaFudge) / n) * 180 / Math.PI;
newLatitude = 0;
for (var i = 0; i < 5 ; ++i) {
newLatitude = (2 * Math.atan(Math.pow(F * a / rho, 1 / n) * Math.pow((1 + e * Math.sin(newLatitude)) / (1 - e * Math.sin(newLatitude)), e / 2))) - Math.PI / 2;
}
newLatitude *= 180 / Math.PI;
return [newLatitude, newLongitude];
}
```

Java:

```
public static Point lambert72toWGS84(double x, double y) {
double newLongitude;
double newLatitude;
double n = 0.77164219;
double F = 1.81329763;
double thetaFudge = 0.00014204;
double e = 0.08199189;
double a = 6378388;
double xDiff = 149910;
double yDiff = 5400150;
double theta0 = 0.07604294;
double xReal = xDiff - x;
double yReal = yDiff - y;
double rho = Math.sqrt(xReal * xReal + yReal * yReal);
double theta = Math.atan(xReal / -yReal);
newLongitude = (theta0 + (theta + thetaFudge) / n) * 180 / Math.PI;
newLatitude = 0;
for (int i = 0; i < 5 ; ++i) {
newLatitude = (2 * Math.atan(Math.pow(F * a / rho, 1 / n) * Math.pow((1 + e * Math.sin(newLatitude)) / (1 - e * Math.sin(newLatitude)), e / 2))) - Math.PI / 2;
}
newLatitude *= 180 / Math.PI;
return new Point(newLatitude, newLongitude);
}
```

C/C++:

```
void LambertToLongLat(double x,double y, double& longt, double& lat){
double n = 0.77164219;
double F = 1.81329763;
double thetaFudge = 0.00014204;
double e= 0.08199189;
double a= 6378388;
double xDiff= 150000;
double yDiff = 5400088.44;
double theta0 = 0.07604294;
double xReal = xDiff-x;
double yReal = yDiff-y;
double rho = sqrt(xReal*xReal + yReal * yReal);
double theta = atan(xReal/-yReal);
longt = (theta0 + (theta+thetaFudge) /n)*180/PI;
lat =0;
for(size_t i=0; i<5; ++i){
lat =(2*atan(pow(F*a/rho,1/n)* pow((1+ e*sin(lat))/(1-e*sin(lat)),e/2)))-PI/2;
}
lat *=180/PI;
}
```

Hope this helps!

Related: if you want to do other transforms (or flip the transorm), check out Transformation between datum’s in JavaScript. Thanks for the tip, Ibrahim (via e-mail)!

**Did this help you out? Like what you see?**

Consider donating.

Consider donating.

I don't run ads on my blog nor do I do this for profit. A donation however would always put a smile on my face though. Thanks!

Hi,

I’ve used your java version for the conversion but unfortunately i am not getting correct latitude and longitude . what is the package that you are using for the Poin class. Any help on this would be appreciated .

Venu

Point is just a basic class with a latitude and a longitude datamember (type: double). That class can easily be created, or you could let the function return an Array instead of a Point instance.

Hi,

Thank you for your email. I’ve applied your algorithm for the below input and i am able to get the output but it isn’t coming as expected as we want to convert Lambert II Étendu Coordinate system (EPSG:27572) to WGS84. Is there any way that we can manipulate the above algorithm that can fit our requirement. Thanks.

Regards,

Venu

Dear Bram,

thanks for this wonderful tool!

I rewrote the conversion in VBA, useful if one has an access database with co-ordinates… The longitudinal works fine, but the lateral is a bit off… Mechelen ends up in Holland…

Am i doing something wrong here?

Cheers,

Wim

Option Compare Database

Function lambert72toWGS84_X(x, y) As Double

Dim newLongitude, newLatitude As Double

Dim n As Double

Dim F As Double

Dim thetaFudge As Double

Dim e As Double

Dim a As Double

Dim xDiff As Double

Dim yDiff As Double

Dim theta0 As Double

Dim xReal As Double

Dim yReal As Double

Dim PI As Double

n = 0.77164219

F = 1.81329763

thetaFudge = 0.00014204

e = 0.08199189

a = 6378388

xDiff = 149910

yDiff = 5400150

theta0 = 0.07604294

xReal = xDiff – x

yReal = yDiff – y

PI = 3.14159265358979

rho = Sqr((xReal * xReal) + (yReal * yReal))

theta = Atn(xReal / -yReal)

newLongitude = (theta0 + ((theta + thetaFudge) / n)) * 180 / PI

lambert72toWGS84_X = newLongitude

End Function

Function lambert72toWGS84_Y(x, y) As Double

Dim newLongitude, newLatitude As Double

Dim n As Double

Dim F As Double

Dim thetaFudge As Double

Dim e As Double

Dim a As Double

Dim xDiff As Double

Dim yDiff As Double

Dim theta0 As Double

Dim xReal As Double

Dim yReal As Double

Dim PI As Double

Dim rho As Double

Dim theta As Double

n = 0.77164219

F = 1.81329763

thetaFudge = 0.00014204

e = 0.08199189

a = 6378388

xDiff = 149910

yDiff = 5400150

theta0 = 0.07604294

PI = 3.141592653

xReal = xDiff – x

yReal = yDiff – y

rho = Sqr((xReal * xReal) + (yReal * yReal))

theta = Atn(xReal / -yReal)

newLatitude = 0

For i = 0 To 4

newLatitude = (2 * (Atn(((F * a) / rho) ^ (1 / n))) * (((1 + (e * Sin(newLatitude))) / (1 – (e * Sin(newLatitude)))) ^ (e / 2))) – (PI / 2)

Next i

newLatitude = newLatitude * (180 / PI)

lambert72toWGS84_Y = newLatitude

End Function

Thanks for sharing, this was really useful for me.

Hello

A first comment about the method for Latitude, you may always sart from an approximate latitude (0,87 for Southern Belgium, 0,89 for Flanders, it reduced the required iteration number and increases accurracy.

Second comment is that your C/C++ script uses slightly different values for xDiff and yDiff

Being involved in studies about management of the territory, source info is aways given in Lambert 72, and points on Google maps are geographical.

I found a page http://www.econet.ulg.ac.be/biogeonet/index.php?pg=951 using a similar script and results are tested and found compliant. But copy-paste method takes a lot of time for maps with more than 500 points to convert.

The formula as proposed by Wim, inserted in an Excel file, give the same result for longitude, difference 1.1 metre being probably due to the reference Google Spherical vs WGS84. But as he points the problem for latitude, I encountered the same, with a difference of 38.65 km, anyway when I tried to write it from Bram’s formulas, and the difference was higher.

Source Lambert 72 141332,132433 (somewhere at the end of Rue des Verreries, Manage). Biogeonet returns 50,502627°N – 004,246559°E.

Excel table returns 51,049445°N – 004,246574° E

Anybody can help ?

New information found on NGI Web site see http://www.ngi.be/Common/Lambert2008/Transformation_Geographic_Lambert_NL.pdf

Values for xDiff and yDiff in C/C++ script is the original, supposed to give co-ordinates as Hayford 1924 ellipsoïd.

Values used in the other scripts are adjusted to match WGS84 co-ordinates. I found better results, matching closely the values of Biogeonet with 149911,1 and 5400149,9

Now about the latitudes

NGI suggest an initial value of NewLatitude.

Calculation for longitude are giving exactly the same results aus your scripts but with another formula not using thetaFudge.

For latitude, an intermediate variable is calculated, based on the initial values.

Using these formulas, I found the expected results.

Thanks for the usefull code Bram.us en Wim.

@Phil De Troy, could you please post the code for the latitude calculation for those of us who are mathematically impaired and not able to translate the formulas into code.

i encounter the same problem…also latitude is way off the expected value.

Does anyone has the code of Phill?

I tried the code of Wim in Access as well and also remarked the problem with the latitude. I wondered if there exists un update of the code. Can I please get it ? I would be very greatful !

Hello. Did anyone found a way to get the latitude working in Excel or Access? The longitude seems to be ok more or less, but the latitude is way off. Thank you very much!

I’ve made an excel with VBA for converting multiple WGS84 coordinates to Lambert 72

download here: https://db.tt/XkcBcYhD

(Thanks to: http://zoologie.umh.ac.be/tc/algorithms.aspx#c2 (Yvan Barbier, UMons, DEMNA))

Axel,

I know it has been a while since you to post the link here. I am at the moment with a problem that I need the same conversion. Would you like to make the excel file back available or send it to mail.

Mvg John

Hi Axel

Your excel file is very useful but I was wondering if there exists a file to convert Lambert 72 into WGS84, the other way around.

Thanks very much.

Kjel Dupon

If you are looking after a tool to convert bunch of coordinates in a CSV file (or in a shapefile), you can always use QGIS. http://qgis.org

You have to import a csv file as a vector layer and specify it is in the Lambert 72 coordinate system (EPSG:31370) and then export it in the WGS84 coordinate system (EPSG:4326)

You are my goddamn hero.

For the delphi/pascal users

Writes the coordinates to a Tpointf (single-precision).

function lambert72toWGS84(const x, y : double): Tpointf;

var XReal, Yreal, rho, theta : double;

i : integer;

const n = 0.77164219;

F = 1.81329763;

thetaFudge = 0.00014204;

e = 0.08199189;

a = 6378388;

xDiff = 149910;

yDiff = 5400150;

theta0 = 0.07604294;

begin

xReal := xDiff – x;

yReal := yDiff – y;

rho := sqrt(xReal * xReal + yReal * yReal);

theta := arctan(xReal/(yReal * -1));

result.X := (theta0 + (theta + thetaFudge)/n) * 180 / system.PI; //Newlongitude

result.Y := 0;

for i := 0 to 4 do

result.Y := (2 * arctan(power(F * a/rho, 1/n) * power((1 + e*sin(result.Y)) / (1-e * sin(result.Y)), e/2))) – System.PI/2;

result.Y := (result.Y * 180)/system.PI; //Newlatitude

end;

Hi,

Thank you for your tutorial,

Please, can you add the code to convert WGS84 to Lambert 72 in java

Thanks