Bienvenue! Inscrivez-vous et rejoignez notre communauté :)
  • Login:

Bienvenue sur Forum SIG - Systèmes d'Information Géographique et Géomatique.

Bienvenue sur le forumSIG. S'il s'agit de votre première visite, assurez vous de faire une recherche préalable dans les FAQ SIG. Vous devez vous inscrire avant de pouvoir poster.

Page 1 sur 4 123 ... DernièreDernière
Affichage des résultats 1 à 15 sur 58
  1. #1

    Date d'inscription
    mai 2006
    Localisation
    Tournai (Belgique)
    Emploi
    Ingénieur d'études
    Organisme
    Cofisun International
    Âge
    35
    Messages
    8

    Question [Données]Formule de conversion coordonnées

    Bonjour,
    je suis nouveau sur le forum donc je m'excuse par avance si le post n'est pas à sa place .


    Je développe une application en Java qui doit recueillir les trames NMEA d'un gps, les analyser pour en tirer les coordonnées, et convertir celles - ci en Lambert II étendu .

    J'ai trouvé presque toutes les informations nécessaires mais il m'en manque encore.
    Voilà comment je procède : je commence par calculer la latitude isométrique, puis je calcule les coordonnées en Lambert II étendu. Toutes les formules que j'utilisent sont tirées du site de l'ign (je suppose donc qu'elles sont justes ). Les voici :

    latitude isométrique : L = ln( tan(PI/4 + phi/2) * ((1-e*sin(phi))/(1+e*sin(phi)))^(e/2) )

    puis : X = Xs + c*exp(-n*L)*sin(n*(lambda - lambda0))
    Y = Ys - c*exp(-n*L)*cos(n*(lambda - lambda0))


    Je connais n, c, Xs, Ys.

    En revanche, pour lambda0 est-ce que je prends le méridien de Paris ou la longitude de Paris par rapport à Greenwich (en gros 0 ou 2°20'14.025... à mettre en radian) ?


    Et deuxième inconnue : e, qui est la première excentricité d'une ellipsoïde de référence... mais laquelle dois-je prendre ? WGS84 ou Clarke 1880 ?

    La logique voudrait que je prenne WGS84, puisque les coordonnées du GPS sont fournies dans ce système...

    L'application que je développe devant être assez précise, j'aimerais éviter de me retrouver à Toulon si je suis à Lille !

    D'avance Merci
    Dernière modification par Phiz ; 19/06/2006 à 14h58.

  2. #2

    Date d'inscription
    février 2005
    Localisation
    Dijon
    Emploi
    Assistant chef de projet - Analyste programmeur
    Organisme
    Atol Conseils & Développements
    Messages
    40

    Par défaut

    un petit coup d'oeil ici :
    http://www.forumsig.org/showthread.php?t=2133

    sinon demande à gily il à l'air de maîtriser le sujet

  3. #3
    Rédacteur honoraire
    Date d'inscription
    octobre 2004
    Localisation
    Trappes
    Emploi
    Ingénieur de recherche
    Organisme
    Laboratoire National de Métrologie et d'Essais
    Âge
    41
    Messages
    186

    Par défaut

    Salut,
    je prends le risque de dire une bêtise, mais pour e je prendrais effecitvement l'excentricité de WGS84. Ce qui me surprend par contre, c'est que tu n'as pas de formule qui fait intervenir les paramètres des deux ellipsoides...
    Et pour lambda0, vue la tête de tes expressions, c'est clairement le méridien de Paris puisque tu veux tes XY en Lambert II.

    J'imagine que quelqu'un pourra confirmer mes dires et lever les doutes!

    JB

  4. #4

    Date d'inscription
    mai 2006
    Localisation
    Tournai (Belgique)
    Emploi
    Ingénieur d'études
    Organisme
    Cofisun International
    Âge
    35
    Messages
    8

    Par défaut

    Salut,

    merci à vous deux.

    Pour les algorithmes, voilà le pdf dans lequel je les ai trouvé :

    http://www.ign.fr/telechargement/MPr...RCE/NTG_71.pdf

    Dans ce document, je me sers des algos 0001 et 0003 ainsi que de la dernière page qui contient les valeurs des constantes (mais ne sont pas toutes là )


    Mais il semblerait quand même que la majorité utilisent ces valeurs pour e et lambda0
    e = 0,08248325676
    lambda0 = 2°20'14.025" = 0,040792344 rad

    A confirmer....

  5. #5

    Date d'inscription
    mai 2006
    Localisation
    Tournai (Belgique)
    Emploi
    Ingénieur d'études
    Organisme
    Cofisun International
    Âge
    35
    Messages
    8

    Lightbulb [Résolu] C'est bon : j'ai la solution !!

    [Post mis à jour le 28 juin 2006]


    Bonjour,

    Juste un dernier petit message pour mettre les choses au point et éclairer tous ceux que çà intéresse.

    Je rappelle le problème : j'ai des coordonnées fournies par un GPS en Degrés Minutes Secondes dans le système WGS84 et je désire les convertir en Lambert 2 étendu (système NTF).

    Voici l'algorithme de conversion que j'ai suivi :
    0) degres-minutes-secondes + orientation (d,m,s,o) radian
    1) coordonnées géographiques WGS84 (phi_w,lambda_w) coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w)
    2) coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w) coordonnées cartésiennes NTF (X_n,Y_n,Z_n)
    3) coordonnées cartésiennes NTF (X_n,Y_n,Z_n) coordonnées géographiques NTF (phi_n,lambda_n)
    4) coordonnées géographiques NTF (phi_n,lambda_n) coordonnées projetées en Lambert II étendu (X_l2e, Y_l2e)
    0)
    Pour la longitude :
    Code:
    lambda_w = d_long + m_long/60 + s_long/3600;
    //Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
    if(orientation == 'W') lambda_w = -1 * lambda_w;
    lambda_w =  lambda_w*Math.PI/180 ;
    Pour la latitude :
    Code:
    phi_w = d_lat + m_lat /60 + s_lat /3600;
    //Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
    if(orientation == 'S') phi_w = -1 * phi_w;
    phi_w =  phi_w*Math.PI/180 ;
    1) J'ai utilisé 2 formules données par l'IGN dans un de leur document ainsi que deux constantes de l'ellipsoide de référence du système WGS84 (les deux demi-axes) :
    Code:
    a_w = 6378137.0;
    b_w = 6356752.314;
    
    //d'où
    e2_w = (a_w*a_w-b_w*b_w)/(a_w*a_w);
    
    //et on a la grande normale de l'ellipsoide WGS84
    N = a_w/sqrt(1-e2_w*pow(sin(phi_w),2)));
    
    //ainsi on obtient
    X_w = N * cos(phi_w) * cos(lambda_w);
    Y_w = N * cos(phi_w) * sin(lambda_w);
    Z_w = N * (1-e2W) * sin(phi_w);
    Ref :
    http://www.ign.fr/telechargement/MPr...RCE/NTG_80.pdf
    http://de.wikipedia.org/wiki/WGS84


    2) Là il n'y a qu'un translation à effectuer :


    Code:
    //on a donc dX = 168.0; dY = 60.0; dZ = -320.0; //et... X_n = X_w + dX; Y_n = Y_w + dY; Z_n = Z_w + dZ;

    Ref :
    http://support.esrifrance.fr/Documents/Generalites/Projections/Generalites/Generalites.htm#2


    3) J'ai utilisé 1 formule donnée par l'IGN toujours dans le même document ainsi que deux constantes de l'ellipsoide de référence du système NTF, Clarke 1880 :
    Code:
    a_n = 6378249.2;
    b_n = 6356515.0;
    
    //d'où
    e2_n = (a_n*a_n-b_n*b_n)/(a_n*a_n);
    
    //on définit une tolérance de convergence
    epsilon = pow(10,-10);
    
    //puis on amorce une boucle de calcul        
    p0=arctan(Z_n/sqrt(X_n*X_n+Y_n*Y_n)*(1-(a_n*e2_n)/(sqrt(X_n*X_n+Y_n*Y_n+Z_n*Z_n))));
    p1=arctan((Z_n/sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*cos(p0))/(sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*pow(sin(p0),2))))));
            
    while(!(Math.abs(p1-p0)<epsilon)){
    p0 = p1; p1 = arctan((Z_n/sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*cos(p0))/(sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*pow(sin(p0),2))))));
    } phi_n = p1; lambda_n = arctan(Y_n/X_n);
    Ref :
    http://www.ign.fr/telechargement/MPr...RCE/NTG_80.pdf
    http://support.esrifrance.fr/Documen...eralites.htm#2


    4) J'utilise les formules de projection et les constantes fournies par l'IGN dans un autre document :
    Code:
    //avant tout on d&#233;finit quelques constantes
    n = 0.7289686274;
    c = 11745793.39;
    Xs = 600000.0;
    Ys = 8199695.768;
        
    e_n = sqrt(e2_n);
    lambda0 = 0.04079234433198;   //correspond &#224; la longitude en radian de Paris (2&#176;20'14.025" E) par rapport &#224; Greenwich
    //puis on calcule la latitude isom&#233;trique
    L = ln(tan(PI/4 + phi_n/2) * pow(((1-e_n*sin(phi_n))/(1+e_n*sin(phi_n))),(e_n/2)));
    
    //enfin on projette
    
    X_l2e = Xs + c*exp((-n*L))*sin(n*(lambdaN-lambda0));
    Y_l2e = Ys - c*.exp((-n*L))*cos(n*(lambdaN-lambda0));
    Ref :
    http://www.ign.fr/telechargement/MPr...RCE/NTG_71.pdf


    J'esp&#232;re que cel&#224; pourra servir &#224; d'autres personnes. En tout cas pour moi, le probl&#232;me est r&#233;solu.
    Dernière modification par Phiz ; 28/06/2006 à 09h12.

  6. #6

    Date d'inscription
    mars 2006
    Emploi
    Responsable d'application
    Âge
    38
    Messages
    61

    Par défaut

    Eh ben, moi ça ne me sert à rien, mais je dois avouer que c'est à la fois agréable et très intéressant d'avoir un suivi aussi précis et détaillé sur ce genre de problème.

    Merci pour tout ceux à qui, comme pour moi, ça permettra d'y voir plus clair.
    InvisiblePinkFlyingSpaghettiMonster prophet.................................................. ..........Pour qu'il y ait le moins de mécontents possible il faut toujours taper sur les mêmes.

  7. #7

    Date d'inscription
    juin 2006
    Messages
    8

    Par défaut

    J'ai pris ton code, j'ai rajouté les 2 3 trucs qui manquaient et j'ai mis tout çà dans une classe Java (vois plus bas)

    Pour l'utiliser, il suffit d'appeler la méthode static comme ceci :

    Code:
    double[]  tabXY = ConversionCoordGPS.WGS84toLambert2e(48,51,39,'N',2,20,43,'E');
    System.out.println("X = " + tabXY[0]);
    System.out.println("Y = " + tabXY[1]);
    Est ce que quelqu'un pourrait l'executer et verifier qu'elle fonctionne bien ?
    Car moi j'ai bien des positions WGS84 à convertir pour tester mais je ne suis pas capable ensuite de vérifier l'exactitude du résultat.

    P'tite question :

    je reçois 47384833,N comme latitude (trame NMEA de mon GPS) et 006506170,E comme longitude.
    Je l'interprete de cette façon sous forme (degré, min, sec) :
    Latitude : 47°38min48sec (orientation N)
    Longitude : 6°50min61sec (orientation E)
    est-correct ? Comment interprêter les 2 derniers chiffres (33 et 70) ?

    Merci d'avance

    Code:
    public class ConversionCoordGPS 
    {
        public final static double[] WGS84toLambert2e(double d_long, double m_long, double s_long, char orientation_long, double d_lat, double m_lat, double s_lat, char orientation_lat)
        {
            double lambda_w, phi_w;
            
            /**************************************************************************************************************/
            /**        0) degres-minutes-secondes + orientation (d,m,s,o) -> radian                                           **/
            /**************************************************************************************************************/
            
            // Pour la longitude
            lambda_w = d_long + m_long/60 + s_long/3600;
            if(orientation_long == 'S') lambda_w = -1 * lambda_w; // Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
            
            // Pour la latitude
            phi_w = d_lat + m_lat /60 + s_lat /3600;
            if(orientation_lat == 'W') phi_w = -1 * phi_w;          // Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
        
            
            /**************************************************************************************************************/
            /**        1) coordonnées géographiques WGS84 (phi_w,lambda_w) -> coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w)  **/
            /**************************************************************************************************************/
            
            // J'ai utilisé 2 formules données par l'IGN dans un de leur document ainsi que deux constantes de 
            // l'ellipsoide de référence du système WGS84 (les deux demi-axes) :
            
            double a_w = 6378137.0;
            double b_w = 6356752.314;
    
            // d'où
            double e2_w = (a_w*a_w-b_w*b_w)/(a_w*a_w);
    
            // et on a la grande normale de l'ellipsoide WGS84
            double N = a_w/Math.sqrt(1 - e2_w * Math.pow(Math.sin(phi_w),2));
    
            // ainsi on obtient
            double X_w = N * Math.cos(phi_w) * Math.cos(lambda_w);
            double Y_w = N * Math.cos(phi_w) * Math.sin(lambda_w);
            double Z_w = N * (1-e2_w) * Math.sin(phi_w);
            
            /**************************************************************************************************************/
            /**        2) coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w) -> coordonnées cartésiennes NTF (X_n,Y_n,Z_n)          **/
            /**************************************************************************************************************/
            
            // Là il n'y a qu'un translation à effectuer :
            
            // on a donc
            double dX = 168.0;
            double dY = 60.0;
            double dZ = -320.0;
    
            // et...
            double X_n = X_w + dX;
            double Y_n = Y_w + dY;
            double Z_n = Z_w + dZ;
            
            /**************************************************************************************************************/
            /**        3) coordonnées cartésiennes NTF (X_n,Y_n,Z_n) -> coordonnées géographiques NTF (phi_n,lambda_n)      **/
            /**************************************************************************************************************/
            
            // J'ai utilisé 1 formule donnée par l'IGN toujours dans le même document ainsi que deux constantes de l'ellipsoide 
            // de référence du système NTF, Clarke 1880 :
            
            double a_n = 6378249.2;
            double b_n = 6356515.0;
    
            // d'où
            double e2_n = (a_n*a_n-b_n*b_n)/(a_n*a_n);
    
            // on définit une tolérance de convergence
            double epsilon = Math.pow(10,-10);
    
            // puis on amorce une boucle de calcul        
            double p0 = Math.atan(Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n)*(1-(a_n*e2_n)/(Math.sqrt(X_n*X_n+Y_n*Y_n+Z_n*Z_n))));
            double p1 = Math.atan((Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*Math.cos(p0))/(Math.sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*Math.pow(Math.sin(p0),2))))));
                    
            while(!(Math.abs(p1-p0)<epsilon))
            {
    
                p0 = p1; p1 = Math.atan((Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*Math.cos(p0))/(Math.sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*Math.pow(Math.sin(p0),2)))))); 
    
            }
            
            double phi_n = p1;
            double lambda_n = Math.atan(Y_n/X_n);
            
            /**********************************************************************************************************************/
            /**        4) coordonnées géographiques NTF (phi_n,lambda_n)  coordonnées projetées en Lambert II étendu (X_l2e, Y_l2e) **/
            /**********************************************************************************************************************/
            
            // J'utilise les formules de projection et les constantes fournies par l'IGN dans un autre document :
            
            // avant tout on définit quelques constantes
            double n = 0.7289686274;
            double c = 11745793.39;
            double Xs = 600000.0;
            double Ys = 8199695.768;
                
            double e_n = Math.sqrt(e2_n);
            double lambda0 = 0.04079234433198;   //correspond à la longitude en radian de Paris (2°20'14.025" E) par rapport à Greenwich
            // puis on calcule la latitude isométrique
            double L = Math.log(Math.tan(Math.PI/4 + phi_n/2) * Math.pow(((1-e_n*Math.sin(phi_n))/(1+e_n*Math.sin(phi_n))),(e_n/2)));
    
            // enfin on projette
    
            double X_l2e = Xs + c*Math.exp((-n*L))*Math.sin(n*(lambda_n-lambda0));
            double Y_l2e = Ys - c*Math.exp((-n*L))*Math.cos(n*(lambda_n-lambda0));
            
            double[] tabXY = new double[2];
            
            tabXY[0] = X_l2e;
            tabXY[1] = Y_l2e;
            
            return tabXY;
        }
    
    }

  8. #8

    Date d'inscription
    mai 2006
    Localisation
    Tournai (Belgique)
    Emploi
    Ingénieur d'études
    Organisme
    Cofisun International
    Âge
    35
    Messages
    8

    Par défaut

    Salut,

    en ce qui concerne les trames NMEA de ton GPS, les donn&#233;es que tu re&#231;ois doivent &#234;tre lues ainsi :
    47384833,N 47&#176; 38,4833' N 47&#176; 38' 28,998" N
    006506170,E 006&#176; 50,6170' E 6&#176; 50' 37,02" E

    Pour plus d'informations sur les trames, tu peux aller voir ici :Pour ce qui est du code que tu as produit, pour le v&#233;rifier, tu peux t&#233;l&#233;charger le logiciel gratuit de l'IGN : Circ&#233;. (http://www.ign.fr/telechargement/MPr...e2000_v133.exe)

    Avec ce logiciel, les coordonn&#233;es que tu as fournies donnent en Lambert II &#233;tendu :
    X = 938479.649
    Y = 2303244.049
    (ATTENTION : L'application des param&#232;tres standard introduit une erreur &#224; ajouter &#224; l'erreur propre sur les coordonn&#233;es initiales (ici de l'ordre de 2 &#224; 5m) ..... c'est pas moi qui le dit, c'est l'IGN .

    En revanche, je viens de m'apercevoir que mon post &#233;tait incomplet... d&#233;sol&#233;. Au tout d&#233;but, lors du passage de degr&#233; vers radian, j'ai oubli&#233; le plus important :

    Code:
    // Pour la longitude
    lambda_w = d_long + m_long/60 + s_long/3600;
    if(orientation_long == 'S') lambda_w = -1 * lambda_w; // Le syst&#232;me de coordonn&#233;es g&#233;ographiques utilis&#233; est postif vers le Nord et vers l'Est
    
     lambda_w =  lambda_w*Math.PI/180 ;
    
            
    // Pour la latitude
    phi_w = d_lat + m_lat /60 + s_lat /3600;
    if(orientation_lat == 'W') phi_w = -1 * phi_w;          // Le syst&#232;me de coordonn&#233;es g&#233;ographiques utilis&#233; est postif vers le Nord et vers l'Est
    
    phi_w = phi_w*Math.PI/180 ;
    Voil&#224;, donc maintenant, ton code est complet, et en le faisant tourner j'obtiens :
    X = 938479.6491473591
    Y = 2303244.050935643
    Le tout &#233;tant donn&#233; en m&#232;tre, tu peux voir par toi m&#234;me que tu as une tr&#232;s l&#233;g&#232;re diff&#233;rence avec ce que donne le logiciel officiel de l'IGN.

    Bon courage pour la suite.
    Il y a deux sortes de diplômés : ceux qui ont appris à apprendre et ceux qui ont appris à penser.

    (Anko Jansen)
    --
    Je suis schizophrène et moi aussi.

    (Thomas Jung)

  9. #9

    Date d'inscription
    juin 2006
    Messages
    8

    Par défaut

    Merci beaucoup !!!! Ca marche nickel maintenant, c'est bien le r&#233;sultat escompt&#233; !

    1) Par contre, ca serait pas plut&#244;t l'inverse pour le test :

    Code:
     
    // Pour la longitude
    if(orientation_long == 'S') -> if(orientation_long == 'W')
    
    // Pour la latitude
    if(orientation_lat == 'W') -> if(orientation_long == 'S')
    La longitude c'est pas E<->W et la latitude N<->S ??
    Enfin ca change RIEN dans mon cas puisque je passe bien les param&#232;tres en latitude -> orientation N et longitude -> orientation E mais pour d'autres personnes ca peut tout changer...
    [B]

    2) Autre chose dont je voudrais &#234;tre s&#251;r :

    La latitude (WGS84) correspond bien au X (Lambert 2 &#233;tendu) et
    la longitude (WGS84) correspond bien au Y (Lambert 2 &#233;tendu) ou c'est l'inverse ?

    J'ai peur de m'&#234;tre un peu m&#233;lang&#233; les pinceaux jusqu'alors

    En tout cas, enore une fois, MERCI 1000 FOIS POUR TON AIDE !!
    Dernière modification par ymerej ; 18/06/2006 à 14h39.

  10. #10

    Date d'inscription
    avril 2005
    Localisation
    France
    Emploi
    Élevé au Pedigree Pal
    Âge
    37
    Messages
    365

    Par défaut

    C'est l'inverse.
    HotShot, péteux d'caillasse en exil.

  11. #11

    Date d'inscription
    mai 2006
    Localisation
    Tournai (Belgique)
    Emploi
    Ingénieur d'études
    Organisme
    Cofisun International
    Âge
    35
    Messages
    8

    Par défaut

    Salut,

    tu as entièrement raison , la longitude c'est bien Est-Ouest et latitude Nord - Sud... j'ai toujours du mal avec les deux notions.. (il m'arrive mais parfois de parler de lontitude !!!! )



    Pour ce qui est des orientations, l'article suivant sur Wikipedia dit bien : "L'axe des X est croissant vers l'Est et l'axe des Y croissant vers le Nord". On peut donc dire qu'à la latitude (qui est bien orientée Nord-Sud), on fait correspondre l'axe des Y, et à la longitude (orientée Est-Ouest) on fait correspondre l'axe des X.


    Pour résumer :
    • Longitude <-> Est-Ouest <-> X
    • Latitude <-> Nord-Sud <-> Y
    Pour finir, une trame NMEA te donne une coordonnée en degré-minute. Pour convertir en degré minute seconde (c'est à dire ° ' "), çà se passe comme avec les heures
    • 1 degré = 60 minutes <-> 1° = 60'
    • 1 minutes = 60 secondes <-> 1' = 60"
    Donc si tu as 47384833 : ca correspond à 47° 38.4833 '. Mais 38.4833 minutes c'est 38 minutes + 0.4833 minutes = 38 minutes + 0.4833*60 secondes = 38 minutes + 28,998 secondes

    d'où le résultat 47384833 = 47° 38' 28.998"

    Dernière remarque pour les trames NMEA, la coordonnée pour la longitude possède 1 chiffre de plus que pour la latitude : celà vient tout simplement du fait que la latitude est forcément un nombre compris entre 0 et 90 (angle mesurer par rapport à l'équateur puis N ou S pour préciser la direction) ; alors que la longitude est un nombre compris entre 0 et 180 (angle mesurer par rapport au Méridien de Greenwich, puis E ou W pour préciser la direction)
    Il y a deux sortes de diplômés : ceux qui ont appris à apprendre et ceux qui ont appris à penser.

    (Anko Jansen)
    --
    Je suis schizophrène et moi aussi.

    (Thomas Jung)

  12. #12

    Date d'inscription
    juin 2006
    Messages
    8

    Par défaut

    Merci Phiz pour toutes tes explications

    Tout est clair désormais.

  13. #13

    Date d'inscription
    juin 2006
    Messages
    8

    Par défaut

    Le code final de la classe est donc :


    Code:
    public class ConversionCoordGPS 
    {
        public final static double[] WGS84toLambert2e(double d_long, double m_long, double s_long, char orientation_long, double d_lat, double m_lat, double s_lat, char orientation_lat)
        {
            double lambda_w, phi_w;
            
            /**************************************************************************************************************/
            /**        0) degres-minutes-secondes + orientation (d,m,s,o) -> radian                                           **/
            /**************************************************************************************************************/
            
            // Pour la longitude
            lambda_w = d_long + m_long/60 + s_long/3600;
            if(orientation_long == 'W') lambda_w = -1 * lambda_w; // Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
    
            lambda_w =  lambda_w*Math.PI/180 ;
    
                    
            // Pour la latitude
            phi_w = d_lat + m_lat /60 + s_lat /3600;
            if(orientation_lat == 'S') phi_w = -1 * phi_w;          // Le système de coordonnées géographiques utilisé est postif vers le Nord et vers l'Est
    
            phi_w = phi_w*Math.PI/180 ;
            
            /**************************************************************************************************************/
            /**        1) coordonnées géographiques WGS84 (phi_w,lambda_w) -> coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w)  **/
            /**************************************************************************************************************/
            
            // J'ai utilisé 2 formules données par l'IGN dans un de leur document ainsi que deux constantes de 
            // l'ellipsoide de référence du système WGS84 (les deux demi-axes) :
            
            double a_w = 6378137.0;
            double b_w = 6356752.314;
    
            // d'où
            double e2_w = (a_w*a_w-b_w*b_w)/(a_w*a_w);
    
            // et on a la grande normale de l'ellipsoide WGS84
            double N = a_w/Math.sqrt(1 - e2_w * Math.pow(Math.sin(phi_w),2));
    
            // ainsi on obtient
            double X_w = N * Math.cos(phi_w) * Math.cos(lambda_w);
            double Y_w = N * Math.cos(phi_w) * Math.sin(lambda_w);
            double Z_w = N * (1-e2_w) * Math.sin(phi_w);
            
            /**************************************************************************************************************/
            /**        2) coordonnées cartésiennes WGS84 (X_w,Y_w,Z_w) -> coordonnées cartésiennes NTF (X_n,Y_n,Z_n)          **/
            /**************************************************************************************************************/
            
            // Là il n'y a qu'un translation à effectuer :
            
            // on a donc
            double dX = 168.0;
            double dY = 60.0;
            double dZ = -320.0;
    
            // et...
            double X_n = X_w + dX;
            double Y_n = Y_w + dY;
            double Z_n = Z_w + dZ;
            
            /**************************************************************************************************************/
            /**        3) coordonnées cartésiennes NTF (X_n,Y_n,Z_n) -> coordonnées géographiques NTF (phi_n,lambda_n)      **/
            /**************************************************************************************************************/
            
            // J'ai utilisé 1 formule donnée par l'IGN toujours dans le même document ainsi que deux constantes de l'ellipsoide 
            // de référence du système NTF, Clarke 1880 :
            
            double a_n = 6378249.2;
            double b_n = 6356515.0;
    
            // d'où
            double e2_n = (a_n*a_n-b_n*b_n)/(a_n*a_n);
    
            // on définit une tolérance de convergence
            double epsilon = Math.pow(10,-10);
    
            // puis on amorce une boucle de calcul        
            double p0 = Math.atan(Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n)*(1-(a_n*e2_n)/(Math.sqrt(X_n*X_n+Y_n*Y_n+Z_n*Z_n))));
            double p1 = Math.atan((Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*Math.cos(p0))/(Math.sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*Math.pow(Math.sin(p0),2))))));
                    
            while(!(Math.abs(p1-p0)<epsilon))
            {
    
                p0 = p1; p1 = Math.atan((Z_n/Math.sqrt(X_n*X_n+Y_n*Y_n))/(1-(a_n*e2_n*Math.cos(p0))/(Math.sqrt((X_n*X_n+Y_n*Y_n)*(1-e2_n*Math.pow(Math.sin(p0),2)))))); 
    
            }
            
            double phi_n = p1;
            double lambda_n = Math.atan(Y_n/X_n);
            
            /**********************************************************************************************************************/
            /**        4) coordonnées géographiques NTF (phi_n,lambda_n)  coordonnées projetées en Lambert II étendu (X_l2e, Y_l2e) **/
            /**********************************************************************************************************************/
            
            // J'utilise les formules de projection et les constantes fournies par l'IGN dans un autre document :
            
            // avant tout on définit quelques constantes
            double n = 0.7289686274;
            double c = 11745793.39;
            double Xs = 600000.0;
            double Ys = 8199695.768;
                
            double e_n = Math.sqrt(e2_n);
            double lambda0 = 0.04079234433198;   //correspond à la longitude en radian de Paris (2°20'14.025" E) par rapport à Greenwich
            // puis on calcule la latitude isométrique
            double L = Math.log(Math.tan(Math.PI/4 + phi_n/2) * Math.pow(((1-e_n*Math.sin(phi_n))/(1+e_n*Math.sin(phi_n))),(e_n/2)));
    
            // enfin on projette
    
            double X_l2e = Xs + c*Math.exp((-n*L))*Math.sin(n*(lambda_n-lambda0));
            double Y_l2e = Ys - c*Math.exp((-n*L))*Math.cos(n*(lambda_n-lambda0));
            
            double[] tabXY = new double[2];
            
            tabXY[0] = X_l2e;
            tabXY[1] = Y_l2e;
            
            return tabXY;
        }
    
    }

  14. #14

    Date d'inscription
    juin 2006
    Localisation
    Cergy-Pontoise
    Emploi
    Etudiant
    Organisme
    E.I.S.T.I.
    Messages
    35

    Par défaut

    Merci, c'est tr&#232;s agr&#233;able d'avoir une bonne explication...

    Je l'ai pass&#233; en C...Si quelqu'un a besoin... (oui, je sais, c'est pas tr&#232;s dur de pass&#233; du java au C, mais bon, au cas ou )

    #edit il faut faire attention &#224; l'utilisation de la fonction "abs" car en C, elle renvoie un entier( ce qui ne nous conviens pas... ). Je vous coneille d'utiliser la fonction "fabs" ;-)
    Dernière modification par Pierre-Alain ; 05/07/2006 à 12h52.

  15. #15

    Date d'inscription
    juillet 2006
    Messages
    2

    Par défaut

    Par le plus grand des hasards, l'un d'entre vous aurait-il implémenté la fonction inverse ?

 

 
Page 1 sur 4 123 ... DernièreDernière

Discussions similaires

  1. [Données] Conversion de coordonnées
    Par adn25fr dans le forum Assistance Technique
    Réponses: 2
    Dernier message: 03/11/2009, 09h32
  2. [Données] Conversion de coordonnées
    Par Flore! dans le forum Assistance Technique
    Réponses: 4
    Dernier message: 27/08/2009, 07h59
  3. [Conversion] Formule de conversion WGS84 > Lambert I
    Par Monsieur Satz dans le forum Assistance Technique
    Réponses: 0
    Dernier message: 17/05/2009, 18h19
  4. [Données] Conversion de coordonnées
    Par Delap dans le forum Assistance Technique
    Réponses: 2
    Dernier message: 15/09/2008, 17h23
  5. [Formule] Coordonnées Lat/Long en X/Y
    Par lebuzz dans le forum Sémiologie et représentation cartographique
    Réponses: 4
    Dernier message: 02/09/2005, 17h03

Liens sociaux

Règles de messages

  • Vous ne pouvez pas créer de nouvelles discussions
  • Vous ne pouvez pas envoyer des réponses
  • Vous ne pouvez pas envoyer des pièces jointes
  • Vous ne pouvez pas modifier vos messages
  •