1

Topic: MySQL Generated Columns OSGB to WGS84

Hi. Does anyone have experience of converting OS grid refs to lat/long using MySQL generated columns functionality? Or is this a foolish thing to attempt? :-)

2

Re: MySQL Generated Columns OSGB to WGS84

Don't know about my SQL, but there are both  SQL Server (in R6 database) and Ms Access function which can do this. Depends though if the Lat/Long you need is WSG84 or OSGB 36.

Mike Weideli

3

Re: MySQL Generated Columns OSGB to WGS84

Hi Mike. I have the same problem lots of people must have had displaying OS grid ref points on Google Maps from a database.

Our OS grid refs are stored in 4 or 6 digit format i.e. SZ0791 or SZ073918.

I know there are client-side JavaScript solutions for point conversion but would prefer to let the database do the processing and store lat/long columns if that's possible.

I believe MySQL only fairly recently introduced 'generated columns' in release 5.7.

4

Re: MySQL Generated Columns OSGB to WGS84

Perhaps someone else can help with this as I don't know a lot about MySQl. R6 holds lat/longs in OSGB36 format.  R6  has functions in Transact SQL which can convert to WSG84 and  there are vb scripts in MS Access available.

Mike Weideli

5

Re: MySQL Generated Columns OSGB to WGS84

Thanks Mike. Assuming they're open source do you have a link to the relevant SQL functions?

6

Re: MySQL Generated Columns OSGB to WGS84

Two UDF below. First converts to WGS84 LL given OS grid refs Eastings and Northings. Second works out Eastings and Northings from a OSGb grid ref. 

USE [NBNData]
GO
/****** Object:  UserDefinedFunction [dbo].[LCOSGTBENtoWGS84LL]    Script Date: 07/17/2017 12:21:49 ******/
SET ANSI_NULLS ON
GO
SET QUOTED_IDENTIFIER ON
GO


/*===========================================================================*\
  Description:    Function to return a WGS84 Lat or Long given a UK E/N
             if not valid returns null
  Parameters:
        @Eastings - UK Eastings in m.
                @Northings - UK Northing in 
                @isLat if 1 returns Lat if 0 then Long
  Returns Lat or Long

  Created:    October 2015

           

  Author: MikeWeideli

\*=========================================================================== */
ALTER FUNCTION [dbo].[LCOSGTBENtoWSG84LL]
  (@Eastings float, @Northings  float, @iSLat bit )

Returns float     
AS
BEGIN
--************************************************************
declare @a float
declare @b float
declare @F0 float
declare @lat0 float
declare @lon0 float
declare @N0 float
declare @E0 float
declare @e2 float
declare @lat float
declare @M float
declare @M1 float
declare @M2 float
declare @M3 float
declare @M4 float
declare @nu float
declare @rho float
declare @eta2 float
declare @seclat float
declare @VII float
declare @VIII float
declare @IX float
declare @X float
declare @XI float
declare @XII float
declare @XIIA float
declare @dE float
declare @lat_1 float
declare @lon_1 float
declare @H float
declare @x_1 float
declare @y_1 float
declare @z_1 float
declare @s float
declare @tx float
declare @ty float
declare @tz float
declare @rxs float
declare @rys float
declare @rzs float
declare @rx float
declare @ry float
declare @rz float
declare @x_2 float
declare @y_2 float
declare @z_2 float
declare @a_2 float
declare @b_2 float
declare @e2_2 float
declare @p float
declare @temp1 float
declare @temp2 float
declare @nu_2 float
declare @LatOld float
declare @Lon float
declare @Long float
declare @Returnvalue float
declare @n float
declare @check float
declare @PI float

set @PI = 3.14159265358979
set @a = 6377563.396
set @b = 6356256.909
set @F0 = 0.9996012717
set @lat0 = 49*@PI/180
set @lon0 = -2*@PI/180
set @N0 = -100000
Set @E0 = 400000
Set @e2 = 1 - (@b*@b)/(@a*@a)
Set @n = (@a-@b)/(@a+@b);
set @lat = @lat0;
set @M = 0;



while (@Northings-@N0-@M >= 0.00001)
Begin

  set @lat = (@Northings-@N0-@M)/(@a*@F0) + @lat
  set @M1 = (1 + @n + (5./4)* power(@n,2) + (5/4)*power(@n,3)) * (@lat-@lat0);
  set @M2 = (3*@n + 3*power(@n,2) + (21/8)*power(@n,3)) * sin(@lat-@lat0) * cos(@lat+@lat0);
  set @M3 = ((15./8)*power(@n,2) + (15/8)*power(@n,3)) * sin(2*(@lat-@lat0)) * cos(2*(@lat+@lat0));
  set @M4 = (35./24)*power(@n,3) * sin(3*(@lat-@lat0)) * cos(3*(@lat+@lat0));
  set @M = @b * @F0 * (@M1 - @M2 + @M3 - @M4) ;

End



Set @nu = @a*@F0/sqrt(1-@e2*power(sin(@lat),2));



set  @rho = @a*@F0*(1-@e2)*power((1-@e2*power(sin(@lat),2)),(-1.5));
set  @eta2 = @nu/@rho-1
set  @secLat = 1./cos(@lat);
set  @VII = tan(@lat)/(2*@rho*@nu);
set  @VIII = tan(@lat)/(24*@rho*power(@nu,3))*(5+3*power(tan(@lat),2)+@eta2-9*power(tan(@lat),2)*@eta2);
set  @IX = tan(@lat)/(720*@rho*power(@nu,5))*(61+90*power(tan(@lat),2)+45*power(tan(@lat),4));
set  @X = @secLat/@nu;
set  @XI = @secLat/(6*power(@nu,3))*(@nu/@rho+2*power(tan(@lat),2));
set  @XII = @secLat/(120*power(@nu,5))*(5+28*power(tan(@lat),2)+24*power(tan(@lat),4));
set  @XIIA = @secLat/(5040*power(@nu,7))*(61+662*power(tan(@lat),2)+1320*power(tan(@lat),4)+720*power(tan(@lat),6));
set  @dE = @Eastings-@E0;

set @lat_1 = @lat - @VII*power(@dE,2) + @VIII*power(@dE,4) - @IX*power(@dE,6);
set @lon_1 = @lon0 + @X*@dE - @XI*power(@dE,3) + @XII*power(@dE,5) - @XIIA*power(@dE,7);



set  @H = 0
set @x_1 = (@nu/@F0 + @H)*cos(@lat_1)*cos(@lon_1);
set @y_1 = (@nu/@F0+ @H)*cos(@lat_1)*sin(@lon_1);
set @z_1 = ((1-@e2)*@nu/@F0 +@H)*sin(@lat_1);




set @s = -20.4894*power(10,-6)
set @tx = 446.448
set @ty = -125.157
set @tz = 542.060
set @rxs = 0.1502
set @rys = 0.2470
set @rzs = 0.8421
set @rx = @rxs*@PI/(180*3600)
set @ry= @rys*@PI/(180*3600)
set @rz = @rzs*@PI/(180*3600)
set @x_2 = @tx + (1+@s)*@x_1 + (-@rz)*@y_1 + (@ry)*@z_1;
set @y_2 = @ty + (@rz)*@x_1 + (1+@s)*@y_1 + (-@rx)*@z_1;
set @z_2 = @tz + (-@ry)*@x_1 + (@rx)*@y_1 + (1+@s)*@z_1;

set @a_2 = 6378137.000
set @b_2 = 6356752.3141
set @e2_2 = 1- (@b_2*@b_2)/(@a_2*@a_2)
set @p = sqrt(power(@x_2,2) + power(@y_2,2));


set @lat = atn2(@z_2,(@p*(1-@e2_2)))
set @latold = 2*@PI;

while (abs(@lat - @latold)>power(10,-16))

Begin
  set @temp1 = @lat
  set @temp2 = @latold
  set @latold = @temp1
  set @lat = @temp2
  set @lat = @latold
  set @latold = @lat
  set @nu_2 = @a_2/sqrt(1-@e2_2*power(sin(@latold),2));
  set @lat = atn2(@z_2+@e2_2*@nu_2*sin(@latold), @p);
End


set @lon = atn2(@y_2,@x_2);
set @H = @p/cos(@lat) - @nu_2;


set @lat = @lat*180/@PI;
set @lon = @lon*180/@PI;


If @isLat = 1
   set @ReturnValue  =  @lat
else
   set @ReturnValue =  @lon

RETURN @ReturnValue

END

USE [NBNData]
GO
/****** Object:  UserDefinedFunction [dbo].[LCOSGBLLToOSGBEN]    Script Date: 07/17/2017 12:27:07 ******/
SET ANSI_NULLS ON
GO
SET QUOTED_IDENTIFIER ON
GO

/*===========================================================================*\
  Description:    Function to return an OSGB Eastings or Northings from a OSGB Lat/Long
             if not valid returns ''
  Parameters:
        @Latitude - latitude in decimal format
                @Longitude in decimal format 
        @Separator text
        @Output 1 = Eastings 2 = Northings 3 = Eastings/Separator/Northings
  Created:    Nov 2014

       

  Author: MikeWeideli

\*=========================================================================== */
ALTER FUNCTION [dbo].[LCOSGBLLToOSGBEN]
(@Latitude float, @Longitude float, @Separator char(1),@Output integer)
RETURNS varchar(20)

AS
BEGIN

Declare @LatW float
Declare @LonW float
Declare @a float
Declare @b float
Declare @f float
Declare @lat0 float
Declare @lon0 float
Declare @n0 float
Declare @e float
Declare @e2 float
Declare @n float
Declare @n2 float
Declare @n3 float
Declare @CosLat float
Declare @SinLat float
Declare @Nu float
Declare @eSq float
Declare @rho float
Declare @eta2 float
Declare @FiveFour float
Declare @FifteenEight float
Declare @TwentyoneEight float
Declare @ThirtyFiveTwentyFour float
Declare @M float
Declare @ma float
Declare @mb float
Declare @mc float
Declare @md float
Declare @lattrueorigin float
Declare @longtrueorigin float

Declare @cos3Lat float
Declare @cos5Lat float
Declare @Tan2Lat float
Declare @Tan4Lat float
Declare @I float
Declare @II float
Declare @III float
Declare @IIIa float
Declare @IV float
Declare @V float
Declare @VI float
Declare @DLon float
Declare @DLat float
Declare @DLon2 float
Declare @DLon3 float
Declare @DLon4 float
Declare @DLon5 float
Declare @DLon6 float
Declare @Northings float
Declare @Eastings float
Declare @Two float
Declare @TwentyFour float
Declare @SevenTwoZero float
Declare @OneTwoZero float
Declare @EastingsInt Integer
Declare @NorthingsInt Integer

Declare @Result varchar(20)

If @Latitude is not null and @Longitude is not null
Begin
    If  @Latitude  = 0 set @Latitude = .0000001
    set  @Lattrueorigin = 49
    set  @Longtrueorigin = -2
    set  @a =  6377563.396
    set  @b =  6356256.908
    set  @f = 0.9996012717
    set  @lat0 = radians(@Lattrueorigin)
    set  @lon0 = radians(@Longtrueorigin)
    set  @latw = radians(@Latitude)
    set  @lonw = radians(@Longitude)


    set  @n0 = -100000
    set  @e = 400000
    set  @e2 = 1 -(@b*@b)/ (@a*@a)
    set  @n = (@a-@b) /(@a+@b)
    set  @n2 = @n*@n
    set  @n3 = @n*@n*@n
    set  @CosLat = cos(@latW)
    set  @SinLat = sin(@latW)
    Set  @Nu = @A*@f / sqrt(1-@e2 * @sinLat * @sinlat)
    set  @rho = @A*@f*(1-@e2) /((1-@e2 * @sinlat) * Power(@sinlat , 1.5))
    set  @eta2 = @Nu/@rho -1
    set  @FiveFour = 1.25
    set  @FifteenEight = 1.875
    set  @TwentyOneEight = 2.625
    Set  @ThirtyFiveTwentyFour = 1.45833333333333
    set  @Two = 2
    set  @TwentyFour = 24
    set  @SevenTwoZero = 720
    set  @OneTwoZero = 120

    set @Ma = (1 + @N + (@FiveFour * @n2) + (@FiveFour * @N3)) * (@latw - @lat0)
    set @Mb = ((3 * @N) + (3 * @N2) + (@TwentyOneEight * @N3)) * Sin(@latw - @lat0) * Cos(@latw + @lat0)
    set @Mc = ((@Fifteeneight * @N2) + (@Fifteeneight * @N3)) * Sin(2 * (@latw - @lat0)) * Cos(2 * (@latw + @lat0))
    set @Md = @ThirtyFiveTwentyFour * @N3 * Sin(3 * (@latw - @lat0)) * Cos(3 * (@latw + @lat0))
    set @M = @b * @F * (@Ma - @Mb + @Mc - @Md)



    set @Cos3lat = @Coslat * @Coslat * @Coslat
    set @Cos5lat = @Cos3lat * @Coslat * @Coslat
    set @tan2lat =  tan(@latw) * tan(@latw)
    set @tan4Lat = @tan2lat * @tan2lat

    set @I = @M + @N0
    set @II = (@Nu / @two) * @SinLat * @Coslat
    set @III = (@Nu / @twentyFour) * @SinLat * @Cos3lat * (5 - @tan2lat + 9 * @Eta2)
    set @IIIA = (@Nu / @SevenTwoZero) * @SinLat * @Cos5lat * (61 - 58 * @tan2lat + @tan4Lat)
    set @IV = @Nu * @Coslat
    set @V = (@Nu / 6) * @Cos3lat * (@Nu / @rho - @tan2lat)
    set @VI = (@Nu / @OneTwoZero) * @Cos5lat * (5 - 18 * @tan2lat + @tan4Lat + 14 * @Eta2 - 58 * @tan2lat * @Eta2)

    set @DLon = @lonw - @lon0
    set @DLon2 = @DLon * @DLon
    set @DLon3 = @DLon2 * @DLon
    set @DLon4 = @DLon3 * @DLon
    set @DLon5 = @DLon4 * @DLon
    set @DLon6 = @DLon5 * @DLon

    set @Northings = @I + @II * @DLon2 + @III * @DLon4 + @IIIA * @DLon6

    set @Eastings =  @E + @IV * @DLon + @V * @DLon3 + @VI * @DLon5 

    set @EastingsInt = ceiling(@Eastings)
    set @NorthingsInt = ceiling(@Northings)


    if @Output = 2
      set @result = cast(@NorthingsInt as  varchar(10))
    else if @output = 1
     set @result = cast(@EastingsInt as varchar(10))
    else
      set @result = cast(@EastingsInt as  varchar(10))  + @Separator +   cast(@NorthingsInt as varchar(10))
end
else
     set @Result = '' 
 
     Return @Result
End

Mike Weideli

7

Re: MySQL Generated Columns OSGB to WGS84

Many thanks Mike. If I find that unworkable with MySQL can you recommend a reliable PHP script as backup?

8

Re: MySQL Generated Columns OSGB to WGS84

Perhaps people are already familiar with this PHP script:

http://www.howtocreate.co.uk/php/gridref.php