// LAMBERT
// http://www.ign.fr/telechargement/MPro/geodesie/CIRCE/NTG_71.pdf

// RGR92 (REUNION)
// http://www.ign.fr/telechargement/MPro/geodesie/CIRCE/systemeReunion.pdf
// http://www.ign.fr/telechargement/MPro/geodesie/CIRCE/NTG_80.pdf

// GAUSS - LABORDE
// http://www.ign.fr/telechargement/MPro/geodesie/CIRCE/NTG_73.pdf

// UTM
// http://www.ign.fr/telechargement/MPro/geodesie/CIRCE/NTG_76.pdf
jQuery.extend({
	
	// constantes UTM : a (demi-grand axe de l'ellipsoïde) et f (n° de fuseau)
	constanteUTMnord: function(a, f){
		return {
			'n'			: 0.9996*a,
			'lambdac'	: 6*f-183,
			'Xs'		: 500000,
			'Ys'		: 0
		};
	},
	constanteUTMsud: function(a, f){
		return {
			'n'			: 0.9996*a,
			'lambdac'	: 6*f-183,
			'Xs'		: 500000,
			'Ys'		: 10000000
		};
	},

	// constantes Lambert France et valeurs associées
	constanteLambert1: function(){
		return {
			'n'		: 0.7604059656,
			'c'		: 11603796.98,
			'Xs'	: 600000,
			'Ys'	: 5657616.674,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378249.2,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},
	constanteLambert2: function(){
		return {
			'n'		: 0.7289686274,
			'c'		: 11745796.39,
			'Xs'	: 600000,
			'Ys'	: 6199695.768,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378249.2,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},
	constanteLambert2Etendu: function(){
		return {
			'n'		: 0.7289686274,
			'c'		: 11745796.39,
			'Xs'	: 600000,
			'Ys'	: 8199695.768,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378249.2,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},
	constanteLambert3: function(){
		return {
			'n'		: 0.6959127966,
			'c'		: 11947992.52,
			'Xs'	: 600000,
			'Ys'	: 6791905.085,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378249.2,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},
	constanteLambert4: function(){
		return {
			'n'		: 0.6712679322,
			'c'		: 12136281.99,
			'Xs'	: 234.358,
			'Ys'	: 7239161.542,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378249.2,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},
	constanteLambert93: function(){
		return {
			'n'		: 0.725607765,
			'c'		: 11754255.426,
			'Xs'	: 700000,
			'Ys'	: 6600000,

			// précision
			'epsilon'	: Math.pow(10,-10),
			// longitude d'origine par rapport au méridien de Greenwich
			'lambdac'	: 0.04079234433,
			// première excentricité de l'ellipsoïde Clarke 1880 français
			'e'			: 0.08248325676,
			// hauteur au-dessus de l'ellipsoïde
			'he'		: 100,
			// demi-grand axe de l'ellipsoïde
			'a'			: 6378137,
			// translation suivant l'axe des x
			'Tx'		: -168,
			// translation suivant l'axe des y
			'Ty'		: -60,
			// translation suivant l'axe des z
			'Tz'		: 320,
			// facteur d'échelle
			'D'			: 0,
			// angles de rotation autour des différents axes, en radians
			'Rx'		: 0,
			'Ry'		: 0,
			'Rz'		: 0
		}
	},

	// Calcul de la latitude isométrique sur un ellipsoïde de première excentricité e au point de latitude phi
	ALGO001: function(phi,e){
		return Math.log( Math.tan ( (Math.PI/4) + (phi/2) ) * Math.pow ((1 - (e * Math.sin(phi))) / (1 + (e * Math.sin(phi))), e/2 ));
	},

	// Calcul de la latitude phi à partir de la lalitude isométrique L
	ALGO002: function(L, e, epsilon){
		var phi	= new Array();
		var i	= 0;
		phi[i]	= 2*Math.atan(Math.exp(L)) - Math.PI/2;
		i++;
		phi[i]	= 2*Math.atan(Math.pow((1+e*Math.sin(phi[i-1]))/(1-e*Math.sin(phi[i-1])), e/2)*Math.exp(L)) - Math.PI/2;

		while(Math.abs(phi[i]-phi[i-1])>=epsilon){
			i++;
			phi[i]	= 2*Math.atan(Math.pow((1+e*Math.sin(phi[i-1]))/(1-e*Math.sin(phi[i-1])), e/2)*Math.exp(L)) - Math.PI/2;
		}
		// phi : latitude
		return phi[i];
	},

	// Transformation de coordonnées géographiques en coordonnées en projection conique conforme de Lambert
	ALGO003: function(lambda, phi, n, c, e, lambdac, Xs, Ys){
		return {
			'X'	: Xs+c*Math.exp(-n*jQuery.ALGO001(phi,e))*Math.sin(n*(lambda-lambdac)),
			'Y'	: Ys-c*Math.exp(-n*jQuery.ALGO001(phi,e))*Math.cos(n*(lambda-lambdac))
		};
	},

	// Transformation de coordonnées en projection conique conforme de Lambert, en coordonnées géographiques.
	ALGO004: function(X, Y, n, c, Xs, Ys, lambdac, e, epsilon){
		var R		= Math.sqrt( Math.pow(X-Xs, 2) + Math.pow(Y-Ys, 2));
		var gamma	= Math.atan((X-Xs)/(Ys-Y));
		var lambda	= lambdac + gamma/n;
		var L		= (-1/n)*Math.log(Math.abs(R/c));

		return {
			'lambda'	: lambda,
			'phi'		: jQuery.ALGO002(L, e, epsilon)
		};
	},

	// transformation de coordonnées géographiques éllipsoïdales en coordonnées cartésiennes
	ALGO009: function(lambda,phi,he,a,e){
		// N : grande normale
		var N	= jQuery.ALGO021(phi, a, e);

		return {
			'X'	: (N + he) * Math.cos(phi) * Math.cos(lambda),
			'Y'	: (N + he) * Math.cos(phi) * Math.sin(lambda),
			'Z'	: (N * (1 - Math.pow(e,2)) + he) * Math.sin(phi)
		};
	},

	// Transformation, pour un ellipsoïde donnée, des coordonnées cartésiennes d'un point en coordonnées géographiques éllipsoïdales
	// par la méthode de Heiskanen-Moritz-Boucher
	ALGO012: function(X,Y,Z,a,e,epsilon){
		var P		= Math.sqrt(Math.pow(X,2)+Math.pow(Y,2));
		var phi		= new Array();
		var i		= 0;
		phi[i]		= Math.atan( Z/ (P * (1 - (a*Math.pow(e,2)) / Math.sqrt(Math.pow(X,2) + Math.pow(Y,2) + Math.pow(Z,2) ) ) ));
		i++;
		var temp	=  Math.pow((1 - ( a * Math.pow(e,2) * Math.cos(phi[i-1])/( P * Math.sqrt(1 - Math.pow(e,2) * Math.pow(Math.sin(phi[i-1]),2))))),-1);
		phi[i]		= Math.atan( temp * Z / P );

		while(Math.abs(phi[i]-phi[i-1]>=epsilon)){
			i++;
			var temp	=  Math.pow((1 - ( a * Math.pow(e,2) * Math.cos(phi[i-1])/( P * Math.sqrt(1 - Math.pow(e,2) * Math.pow(Math.sin(phi[i-1]),2))))),-1);
			phi[i]		= Math.atan( temp * Z / P );
		}

		return {
			'lambda'	: Math.atan(Y/X),
			'phi'		: phi[i],
			'he'		: (P/Math.cos(phi[i])) - (a/Math.sqrt(1 - Math.pow(e,2) * Math.pow(Math.sin(phi[i]),2)))
		};
	},

	// A partir d'un jeu de 7 paramètres (3 translations, 1 facteur d'échelle et 3 rotations) de passage du système (1)
	// vers le système (2), et des coordonnées cartésiennes tridimensionnelles dans le système (1), calcul des coordonnées cartésiennes
	// tridimenseionnelles dans le système (2)
	ALGO013: function(Tx,Ty,Tz,D,Rx,Ry,Rz,X,Y,Z){
		return {
			'X'	: Tx + X * (1 + D) + Z * Ry - Y * Rz,
			'Y'	: Ty + Y * (1 + D) + X * Rz - Z * Rx,
			'Z'	: Tz + Z * (1 + D) + Y * Rx - X * Ry
		};
	},

	// Calcul de la grande normale de l'ellipsoïde
	ALGO021: function(phi, a, e){
		return a/Math.sqrt(1-Math.pow(e,2)*Math.pow(Math.sin(phi),2) );
	},

	// conversion de lambert vers WGS-84
	lambert_WGS84: function(systeme, X, Y){
		var c	= {};
		switch(systeme){
			case "LambertI":
				c	= jQuery.constanteLambert1();
				break;
			case "LambertII":
				c	= jQuery.constanteLambert2();
				break;
			case "LambertIIE":
				c	= jQuery.constanteLambert2Etendu();
				break;
			case "LambertIII":
				c	= jQuery.constanteLambert3();
				break;
			case "LambertIV":
				c	= jQuery.constanteLambert4();
				break;
			case "Lambert93":
				c	= jQuery.constanteLambert93();
				break;
		}
		var coords	= jQuery.ALGO004(X, Y, c.n, c.c, c.Xs, c.Ys, c.lambdac, c.e, c.epsilon);

		coords		= jQuery.ALGO009(coords.lambda,coords.phi,c.he,c.a,c.e);

		coords		= jQuery.ALGO013(c.Tx,c.Ty,c.Tz,c.D,c.Rx,c.Ry,c.Rz,coords.X,coords.Y,coords.Z);

		// ellipsoïdes WGS84
		var a		= 6378137;
		var f		= 1/298.257223563;
		var b		= a*(1-f);
		var e		= Math.sqrt((Math.pow(a,2) - Math.pow(b,2))/Math.pow(a,2));

		coords = jQuery.ALGO012(coords.X,coords.Y,coords.Z,a,e,c.epsilon);

		return {
			'lambda'	: jQuery.rad2deg(coords.lambda),
			'phi'		: jQuery.rad2deg(coords.phi)
		};
	},

	// conversion de WGS-84 vers lambert
	WGS84_lambert: function(systeme, lambda, phi){
		lambda		= jQuery.deg2rad(lambda);
		phi			= jQuery.deg2rad(phi);

		// ellipsoïdes WGS84
		var c		= {};
		c.a			= 6378137;
		var f		= 1/298.257223563;
		var b		= c.a*(1-f);
		c.e			= Math.sqrt((Math.pow(c.a,2) - Math.pow(b,2))/Math.pow(c.a,2));
		c.he		= 0;

		var coords	= jQuery.ALGO009(lambda,phi,c.he,c.a,c.e);

		coords.X	+= 168;
		coords.Y	+= 60;
		coords.Z	+= -320;

		switch(systeme){
			case "LambertI":
				c	= jQuery.constanteLambert1();
				break;
			case "LambertII":
				c	= jQuery.constanteLambert2();
				break;
			case "LambertIIE":
				c	= jQuery.constanteLambert2Etendu();
				break;
			case "LambertIII":
				c	= jQuery.constanteLambert3();
				break;
			case "LambertIV":
				c	= jQuery.constanteLambert4();
				break;
			case "Lambert93":
				c	= jQuery.constanteLambert93();
				break;
		}

		coords	= jQuery.ALGO012(coords.X,coords.Y,coords.Z,c.a,c.e,c.epsilon);

		coords	= jQuery.ALGO003(coords.lambda, coords.phi, c.n, c.c, c.e, c.lambdac, c.Xs, c.Ys);

		return {
			'X'	: coords.X,
			'Y'	: coords.Y
		};
	},

	// converti des radians en degrés
	rad2deg: function(radians){
		return 360 * radians/(2 * Math.PI);
	},

	// converti des degrés en radians
	deg2rad: function(degrees){
		return (2 * Math.PI * degrees)/360;
	}
});