Deepseek responses when asked to generate a function to convert OSGB36 to WGS84 co-ordinates using pure Ruby

01.prompt.md

Generate a Ruby function to convert UK Ordnance Survey OSGB36 co-ordinates to WGS84 latitude and longitude just using pure Ruby and not using a gem.

02.response.md

To convert UK Ordnance Survey OSGB36 coordinates (easting and northing) to WGS84 latitude and longitude, you can use a pure Ruby implementation of the coordinate transformation. This involves several steps, including converting OSGB36 easting/northing to latitude/longitude in the OSGB36 datum, and then transforming those coordinates to WGS84.

Here’s a Ruby function to perform this conversion:

def osgb36_to_wgs84(easting, northing)
  # OSGB36 parameters
  a = 6377563.396
  b = 6356256.909
  e0 = 400000.0
  n0 = -100000.0
  f0 = 0.9996012717
  phi0 = 49.0 * Math::PI / 180.0
  lambda0 = -2.0 * Math::PI / 180.0

  # WGS84 parameters
  a_wgs84 = 6378137.0
  b_wgs84 = 6356752.314245
  e2_wgs84 = (a_wgs84**2 - b_wgs84**2) / a_wgs84**2

  # Airy 1830 ellipsoid parameters
  e2 = (a**2 - b**2) / a**2
  n = (a - b) / (a + b)

  # Convert easting and northing to OSGB36 latitude and longitude
  phi_prime = ((northing - n0) / (a * f0)) + phi0
  m = calculate_m(phi_prime, a, e2)

  while (northing - n0 - m).abs >= 0.00001
    phi_prime += (northing - n0 - m) / (a * f0)
    m = calculate_m(phi_prime, a, e2)
  end

  nu = a * f0 / Math.sqrt(1 - e2 * Math.sin(phi_prime)**2)
  rho = a * f0 * (1 - e2) / (1 - e2 * Math.sin(phi_prime)**2)**1.5
  eta2 = nu / rho - 1

  tan_phi = Math.tan(phi_prime)
  tan2_phi = tan_phi**2
  tan4_phi = tan2_phi**2

  sec_phi = 1.0 / Math.cos(phi_prime)
  nu3 = nu**3
  nu5 = nu**5

  vii = tan_phi / (2 * rho * nu)
  viii = tan_phi / (24 * rho * nu3) * (5 + 3 * tan2_phi + eta2 - 9 * tan2_phi * eta2)
  ix = tan_phi / (720 * rho * nu5) * (61 + 90 * tan2_phi + 45 * tan4_phi)

  x = (easting - e0) / nu
  x2 = x**2
  x3 = x**3
  x4 = x2**2
  x5 = x3 * x2

  phi = phi_prime - vii * x2 + viii * x4 - ix * x5
  lambda = lambda0 + (x / Math.cos(phi_prime)) - (x3 / (6 * nu2 * Math.cos(phi_prime))) * (1 + 2 * tan2_phi + eta2) + (x5 / (120 * nu4 * Math.cos(phi_prime))) * (5 + 28 * tan2_phi + 24 * tan4_phi)

  # Convert OSGB36 latitude and longitude to WGS84
  dx = -446.448
  dy = 125.157
  dz = -542.060
  s = 20.4894 / 1_000_000.0

  # Helmert transformation
  lambda_rad = lambda
  phi_rad = phi

  sin_phi = Math.sin(phi_rad)
  cos_phi = Math.cos(phi_rad)
  sin_lambda = Math.sin(lambda_rad)
  cos_lambda = Math.cos(lambda_rad)

  x_osgb36 = (nu + h) * cos_phi * cos_lambda
  y_osgb36 = (nu + h) * cos_phi * sin_lambda
  z_osgb36 = ((1 - e2) * nu + h) * sin_phi

  x_wgs84 = x_osgb36 + dx
  y_wgs84 = y_osgb36 + dy
  z_wgs84 = z_osgb36 + dz

  # Convert Cartesian to WGS84 latitude and longitude
  p = Math.sqrt(x_wgs84**2 + y_wgs84**2)
  lambda_wgs84 = Math.atan2(y_wgs84, x_wgs84)
  phi_wgs84 = Math.atan2(z_wgs84, p * (1 - e2_wgs84))

  # Convert radians to degrees
  latitude_wgs84 = phi_wgs84 * 180.0 / Math::PI
  longitude_wgs84 = lambda_wgs84 * 180.0 / Math::PI

  [latitude_wgs84, longitude_wgs84]
end

def calculate_m(phi, a, e2)
  a * ((1 - e2 / 4 - 3 * e2**2 / 64 - 5 * e2**3 / 256) * phi -
       (3 * e2 / 8 + 3 * e2**2 / 32 + 45 * e2**3 / 1024) * Math.sin(2 * phi) +
       (15 * e2**2 / 256 + 45 * e2**3 / 1024) * Math.sin(4 * phi) -
       (35 * e2**3 / 3072) * Math.sin(6 * phi))
end

# Example usage:
easting = 651409.903
northing = 313177.270
latitude, longitude = osgb36_to_wgs84(easting, northing)
puts "Latitude: #{latitude}, Longitude: #{longitude}"

Explanation:

1. OSGB36 to Latitude/Longitude Conversion: The function first converts the OSGB36 easting and northing to latitude and longitude in the OSGB36 datum using the Transverse Mercator projection.
2. Helmert Transformation: The OSGB36 latitude and longitude are then transformed to WGS84 using a Helmert transformation, which accounts for the differences in the ellipsoids and datums.
3. Cartesian to WGS84 Conversion: Finally, the Cartesian coordinates are converted back to WGS84 latitude and longitude.

Notes:

• The Helmert transformation parameters (dx, dy, dz, s) are specific to the OSGB36 to WGS84 transformation.
• The function assumes the height (h) is zero. If you have height data, you can include it in the calculation.

This implementation is a simplified version and may not be as accurate as using a specialized library, but it should work for basic conversions.