A ruby class extracted from an internal lib that transforms UTM coords to Geodesic and from Asia coords to UTM. It's specific to deal with The City Council of Barcelona's opened data.

gistfile1.rb

# encoding: utf-8
##
# Canonizes a geoposition as a geodesic coordinate.
#
# The “Asia” coords were created by the City Council of Barcelona
# and they are still used in some internal systems. The problem
# appears when http://opendata.bcn.cat/ opens data using this non
# standard system.
#
# @author Arnau Siches
class Coords
  # @return [Array] the pair +[latitude, longitude]+.
  attr_reader :value
  attr_reader :latitude, :longitude

  # @return [Hash] the original value (+:value+) Array and the original type
  #   (+:type+).
  attr_reader :original

  ##
  # @param [Array] coords The pair (x, y).
  #
  # @param [Symbol] type The type of the coords (:asia, :utm, :geo)
  def initialize(coords, type)
    @original = {
      :value => coords,
      :type => type
    }
    @type = type

    @value = to_geo(coords)
    @latitude = @value[0]
    @longitude = @value[1]
  end

  private

  ##
  # Converts Asia coords to UTM coords.
  #
  # @param [Array] The coords in Asia.
  #
  # @return [Array] The coords in UTM.
  def asia_to_utm(coords)
    x = (coords[0] + 400000000.0) / 1000
    y = (coords[1] + 4500000000.0) / 1000
    @type = :utm
    return [x, y]
  end

  ##
  # Converts UTM or Asia coords to Geodesic coords.
  #
  # @param [Array] The coords in UTM or Asia.
  #
  # @return [Array] The coords in Geodesic.
  def to_geo(coords)
    case @type
    when :asia
      to_geo(asia_to_utm(coords))
    when :utm
      utm_to_geo(coords)
    end
  end

  ##
  # Converts UTM coords to Geodesic coords.
  #
  # @param [Array] The coords in UTM.
  #
  # @return [Array] The coords in Geodesic.
  def utm_to_geo(coords)
    x = coords[0]
    y = coords[1]

    tx = -131.250417839
    ty = -206.9696618318
    exy = (1.3219081270/1000000)
    ang = (-1.6446198/3600) * Math::PI / 180

    # Coords to UTM ETRS89
    tmpx = tx + ((1 + exy) * ((x * Math.cos(ang)) - (y * Math.sin(ang))))
    tmpy = ty + ((1 + exy) * ((x * Math.sin(ang)) + (y * Math.cos(ang))))

    major_semiaxe = 6378137.0
    minor_semiaxe = 6356752.31414035
    second_excentricity = ((major_semiaxe**2 - minor_semiaxe**2)**0.5) / (minor_semiaxe)
    second_excentricity_power_2 = second_excentricity**2
    curvature_polar_radius = major_semiaxe**2 / minor_semiaxe

    x_recoil = tmpx - 500000.0

    decimal_longitude = (31 * 6) - 183.0
    decimal_latitude = tmpy / (6366197.724 * 0.9996)

    nu = curvature_polar_radius * 0.9996 / (1 + second_excentricity_power_2 * Math.cos(decimal_latitude)**2)**0.5
    a = x_recoil / nu

    a1 = Math.sin(2 * decimal_latitude)
    a2 = a1 * Math.cos(decimal_latitude)**2
    j2 = decimal_latitude + (a1 / 2)
    j4 = ((3 * j2) + a2) / 4
    j6 = ((5 * j4) + a2 * Math.cos(decimal_latitude)**2) / 3

    alfa = 3 * second_excentricity_power_2 / 4
    beta = 5 * alfa**2 / 3
    gamma = 35 * alfa**3 / 27
    bsubzeta = 0.9996 * curvature_polar_radius * (decimal_latitude - alfa * j2 + beta * j4 - gamma * j6)

    b = (tmpy - bsubzeta) / nu
    dseta = (second_excentricity_power_2 * a**2 * Math.cos(decimal_latitude)**2) / 2
    xi = a * (1 - (dseta / 3))
    eta = b * (1 - dseta) + decimal_latitude

    ee = Math.exp(1.0)
    sin_h_xi = (ee**xi - ee**-xi) / 2
    delta_lambda = Math.atan(sin_h_xi / Math.cos(eta))
    tau = Math.atan(Math.cos(delta_lambda) * Math.tan(eta))

    longitude = ((delta_lambda / Math::PI) * 180) + decimal_longitude
    latitude = ((decimal_latitude + (((1 + (second_excentricity_power_2 * Math.cos(decimal_latitude)**2) - (3 / 2 * second_excentricity_power_2 * Math.sin(decimal_latitude) * Math.cos(decimal_latitude) * (tau - decimal_latitude)))) * (tau - decimal_latitude))) / Math::PI) * 180

    @type = :geo

    return [latitude, longitude]
  end
end