var deg2rad = Math.PI / 180;
var rad2deg = 180.0 / Math.PI;
var pi = Math.PI;
var deciLat;
var deciLon;


// Processes WGS84 lat and lon in decimal form with S and W as -ve
function Transform(WGlat, WGlon, height) {
	//first off convert to radians
	var radWGlat = WGlat * deg2rad;
	var radWGlon = WGlon * deg2rad;
	
	//these are the values for WGS86(GRS80) to OSGB36(Airy)
	var a = 6378137;              // WGS84_AXIS
	var e = 0.00669438037928458;  // WGS84_ECCENTRIC
	var h = height;               // height above datum  (from GPS GGA sentence)
	var a2 = 6377563.396;         //OSGB_AXIS
	var e2 = 0.0066705397616;     // OSGB_ECCENTRIC 
	var xp = -446.448;
	var yp = 125.157;
	var zp = -542.06;
	var xr = -0.1502;
	var yr = -0.247;
	var zr = -0.8421;
	var s = 20.4894;

	// convert to cartesian; lat, lon are radians
	var sf = s * 0.000001;
	var v = a / (Math.sqrt(1 - (e * (Math.sin(radWGlat) * Math.sin(radWGlat)))));
	var x = (v + h) * Math.cos(radWGlat) * Math.cos(radWGlon);
	var y = (v + h) * Math.cos(radWGlat) * Math.sin(radWGlon);
	var z = ((1 - e) * v + h) * Math.sin(radWGlat);

	// transform cartesian
	var xrot = (xr / 3600) * deg2rad;
	var yrot = (yr / 3600) * deg2rad;
	var zrot = (zr / 3600) * deg2rad;
	var hx = x + (x * sf) - (y * zrot) + (z * yrot) + xp;
	var hy = (x * zrot) + y + (y * sf) - (z * xrot) + yp;
	var hz = (-1 * x * yrot) + (y * xrot) + z + (z * sf) + zp;

	// Convert back to lat, lon
	var newLon = Math.atan(hy / hx);
	var p = Math.sqrt((hx * hx) + (hy * hy));
	var newLat = Math.atan(hz / (p * (1 - e2)));
	v = a2 / (Math.sqrt(1 - e2 * (Math.sin(newLat) * Math.sin(newLat))));
	var errvalue = 1.0;
	var lat0 = 0;
	while (errvalue > 0.001) {
		lat0 = Math.atan((hz + e2 * v * Math.sin(newLat)) / p);
		errvalue = Math.abs(lat0 - newLat);
		newLat = lat0;
	}

	//convert back to degrees
	newLat = newLat * rad2deg;
	newLon = newLon * rad2deg;
	
	//convert lat and lon (OSGB36)  to OS 6 figure northing and easting
	return LLtoNE(newLat, newLon);
}

//converts lat and lon (OSGB36)  to OS 6 figure northing and easting
function LLtoNE(lat, lon) {
	var phi = lat * deg2rad;          // convert latitude to radians
	var lam = lon * deg2rad;          // convert longitude to radians
	var a = 6377563.396;              // OSGB semi-major axis
	var b = 6356256.91;               // OSGB semi-minor axis
	var e0 = 400000;                  // easting of false origin
	var n0 = -100000;                 // northing of false origin
	var f0 = 0.9996012717;            // OSGB scale factor on central meridian
	var e2 = 0.0066705397616;         // OSGB eccentricity squared
	var lam0 = -0.034906585039886591; // OSGB false east
	var phi0 = 0.85521133347722145;   // OSGB false north
	var af0 = a * f0;
	var bf0 = b * f0;

	// easting
	var slat2 = Math.sin(phi) * Math.sin(phi);
	var nu = af0 / (Math.sqrt(1 - (e2 * (slat2))));
	var rho = (nu * (1 - e2)) / (1 - (e2 * slat2));
	var eta2 = (nu / rho) - 1;
	var p = lam - lam0;
	var IV = nu * Math.cos(phi);
	var clat3 = Math.pow(Math.cos(phi), 3);
	var tlat2 = Math.tan(phi) * Math.tan(phi);
	var V = (nu / 6) * clat3 * ((nu / rho) - tlat2);
	var clat5 = Math.pow(Math.cos(phi), 5);
	var tlat4 = Math.pow(Math.tan(phi), 4);
	var VI = (nu / 120) * clat5 * ((5 - (18 * tlat2)) + tlat4 + (14 * eta2) - (58 * tlat2 * eta2));
	var east = e0 + (p * IV) + (Math.pow(p, 3) * V) + (Math.pow(p, 5) * VI);

	// northing
	var n = (af0 - bf0) / (af0 + bf0);
	var M = Marc(bf0, n, phi0, phi);
	var I = M + (n0);
	var II = (nu / 2) * Math.sin(phi) * Math.cos(phi);
	var III = ((nu / 24) * Math.sin(phi) * Math.pow(Math.cos(phi), 3)) * (5 - Math.pow(Math.tan(phi), 2) + (9 * eta2));
	var IIIA = ((nu / 720) * Math.sin(phi) * clat5) * (61 - (58 * tlat2) + tlat4);
	var north = I + ((p * p) * II) + (Math.pow(p, 4) * III) + (Math.pow(p, 6) * IIIA);

	// make whole number values
	east = Math.round(east);   // round to whole number
	north = Math.round(north); // round to whole number

	return east + "," + north;
}

// a function used in LLtoNE  - that's all I know about it
function Marc(bf0, n, phi0, phi) {
	return bf0 * (((1 + n + ((5 / 4) * (n * n)) + ((5 / 4) * (n * n * n))) * (phi - phi0))
	   - (((3 * n) + (3 * (n * n)) + ((21 / 8) * (n * n * n))) * (Math.sin(phi - phi0)) * (Math.cos(phi + phi0)))
	   + ((((15 / 8) * (n * n)) + ((15 / 8) * (n * n * n))) * (Math.sin(2 * (phi - phi0))) * (Math.cos(2 * (phi + phi0))))
	   - (((35 / 24) * (n * n * n)) * (Math.sin(3 * (phi - phi0))) * (Math.cos(3 * (phi + phi0)))));
}

function transform(lat, lon, a, e, h, a2, e2, xp, yp, zp, xr, yr, zr, s)
{
  // convert to cartesian; lat, lon are radians
  sf = s * 0.000001;
  v = a / (sqrt(1 - (e *(sin(lat) * sin(lat)))));
  x = (v + h) * cos(lat) * cos(lon);
  y = (v + h) * cos(lat) * sin(lon);
  z = ((1 - e) * v + h) * sin(lat);

  xrot = (xr / 3600) * deg2rad;
  yrot = (yr / 3600) * deg2rad;
  zrot = (zr / 3600) * deg2rad;
  
  hx = x + (x * sf) - (y * zrot) + (z * yrot) + xp;
  hy = (x * zrot) + y + (y * sf) - (z * xrot) + yp;
  hz = (-1 * x * yrot) + (y * xrot) + z + (z * sf) + zp;

  // Convert back to lat, lon
  lon = atan(hy / hx);
  p = sqrt((hx * hx) + (hy * hy));
  lat = atan(hz / (p * (1 - e2)));
  v = a2 / (sqrt(1 - e2 * (sin(lat) * sin(lat))));
  errvalue = 1.0;
  lat0 = 0;
  while (errvalue > 0.001)
  {
    lat0 = atan((hz + e2 * v * sin(lat)) / p);  
    errvalue = abs(lat0 - lat);
    lat = lat0;
  }
  h = p / cos(lat) - v;
  var geo = { latitude: lat, longitude: lon };
  return(geo);
}
//===================================================================
function ne2ll(east , north)  // (east, north, grid)
{
// converts NGR easting and nothing to lat, lon.
// input metres, output radians
  var grid = 1;
  var nX = Number(north);
  var eX = Number(east);
  var locres = 0;

  // validate - check all are numbers and something is entered
  if ((String(nX) == "NaN") || (String(eX) == "NaN") || (north.length == 0) || (east.length == 0))
  {
    alert("Invalid northings or eastings");
    return;
  }  
     
    a = 6377563.396;       // OSI semi-major
    b = 6356256.91;        // OSI semi-minor
    e0 = 400000;           // easting of false origin
    n0 = -100000;          // northing of false origin
    f0 = 0.9996012717;     // OSI scale factor on central meridian
    e2 = 0.0066705397616;  // OSI eccentricity squared
    lam0 = -0.034906585039886591;  // OSI false east
    phi0 = 0.85521133347722145;    // OSI false north
    //en2ngr(round(eX), round(nX), "British");
 
  var af0 = a * f0;
  var bf0 = b * f0;
  var n = (af0 - bf0) / (af0 + bf0);
  var Et = east - e0;
  var phid = InitialLat(north, n0, af0, phi0, n, bf0);
  
  var nu = af0 / (sqrt(1 - (e2 * (sin(phid) * sin(phid)))));
  var rho = (nu * (1 - e2)) / (1 - (e2 * (sin(phid)) * (sin(phid))));
  var eta2 = (nu / rho) - 1;
  var tlat2 = tan(phid) * tan(phid);
  var tlat4 = pow(tan(phid), 4);
  var tlat6 = pow(tan(phid), 6);
  var clatm1 = pow(cos(phid), -1);
  var VII = tan(phid) / (2 * rho * nu);
  var VIII = (tan(phid) / (24 * rho * (nu * nu * nu))) * (5 + (3 * tlat2) + eta2 - (9 * eta2 * tlat2));
  var IX = ((tan(phid)) / (720 * rho * pow(nu, 5))) * (61 + (90 * tlat2) + (45 * tlat4));
  var phip = (phid - ((Et * Et) * VII) + (pow(Et, 4) * VIII) - (pow(Et, 6) * IX)); 
  var X = pow(cos(phid), -1) / nu;
  var XI = (clatm1 / (6 * (nu * nu * nu))) * ((nu / rho) + (2 * (tlat2)));
  var XII = (clatm1 / (120 * pow(nu, 5))) * (5 + (28 * tlat2) + (24 * tlat4));
  var XIIA = clatm1 / (5040 * pow(nu, 7)) * (61 + (662 * tlat2) + (1320 * tlat4) + (720 * tlat6));
  var lambdap = (lam0 + (Et * X) - ((Et * Et * Et) * XI) + (pow(Et, 5) * XII) - (pow(Et, 7) * XIIA));


  var geo = convert_to_wgs(grid, phip, lambdap);
  
  var lat = geo.latitude * rad2deg;
  var lon = geo.longitude * rad2deg;
  return new GLatLng(lat,  lon);
  //alert(lat+','+lon);
}  
//=================================================================== 
function InitialLat(north, n0, af0, phi0, n, bf0)
{
  var phi1 = ((north - n0) / af0) + phi0;
  var M = Marc(bf0, n, phi0, phi1);
  var phi2 = ((north - n0 - M) / af0) + phi1;
  var ind = 0;
  while ((abs(north - n0 - M) > 0.00001) && (ind < 20))  // max 20 iterations in case of error
  {  
	ind = ind + 1;
	phi2 = ((north - n0 - M) / af0) + phi1;
    M = Marc(bf0, n, phi0, phi2);
    phi1 = phi2;
  }
  return(phi2);  
}
//===================================================================
function convert_to_wgs(grid, phip, lambdap)
{
  var WGS84_AXIS = 6378137;
  var WGS84_ECCENTRIC = 0.00669438037928458;
  var OSGB_AXIS = 6377563.396;
  var OSGB_ECCENTRIC = 0.0066705397616;
  var IRISH_AXIS = 6377340.189;
  var IRISH_ECCENTRIC = 0.00667054015;
  var INT24_AXIS = 6378388.000;
  var INT24_ECCENTRIC = 0.0067226700223333;
  var height = 10;  // dummy height
  if (grid == 1)
  {
    var geo = transform(phip, lambdap, OSGB_AXIS, OSGB_ECCENTRIC, height, WGS84_AXIS, WGS84_ECCENTRIC, 446.448, -125.157, 542.06, 0.1502, 0.247, 0.8421, 20.4894);
  }
  if (grid == 2)
  {
    var geo = transform(phip, lambdap, IRISH_AXIS, IRISH_ECCENTRIC, height, WGS84_AXIS, WGS84_ECCENTRIC, 482.53, -130.596, 564.557, -1.042, -0.214, -0.631, -8.15);
  }
  if (grid == 3)
  {  
    var geo = transform(phip, lambdap, INT24_AXIS, INT24_ECCENTRIC, height,  WGS84_AXIS, WGS84_ECCENTRIC, -83.901, -98.127, -118.635, 0, 0, 0, 0);
  } 
  return(geo);
}
//-------------------------------------------------------------------
// Maths functions

function mod(y, x)
{
if (y >= 0)
  return y - x * floor(y / x);
else
  return y + x * (floor(-y / x) + 1.0);
}

function atan2(y, x)
{
	return Math.atan2(y, x);
}

function sqrt(x)
{
	return Math.sqrt(x);
}

function tan(x)
{
	return Math.tan(x);
}

function sin(x)
{
	return Math.sin(x);
}

function cos(x)
{
	return Math.cos(x);
}

function acos(x)
{
	return Math.acos(x);
}

function floor(x)
{
	return Math.floor(x);
}

function round(x)
{
	return Math.round(x);
}

function ceil(x)
{
	return Math.ceil(x)
}

function ln(x)
{
	return Math.log(x);
}

function abs(x)
{
	return Math.abs(x);
}

function pow(x, y)
{
	return Math.pow(x, y);
}

function atan(x)
{
	return Math.atan(x);
}

function chr(x)
{
return String.fromCharCode(x);
}

function round(x)
{
	return Math.round(x);
}