dimanche 16 juin 2013

Distance entre deux codes postaux au Canada et radius

Code postal

Le code postal au Canada permet de trouver l'emplacement d'un bâtiment d'une façon assez précise.
Il y en a près de 1 million. 
Une base de données des codes postaux Canadien peut être téléchargé à cette adresse: http://geocoder.ca/?freedata=1

Le fichier téléchargé est un fichier csv, il a la structure suivante

Code postal Latitude Longitude Ville Province
J5R5S8 45.39996 -73.506515 CandiacQC
J5R5S9 45.422333 -73.488784 La Prairie QC


Chaque code postal possède une longitude et une latitude. Ces informations peut être utile afin de connaitre la distance d'un lieu par rapport à un autre. Par exemple après avoir retourné tous les items d'une recherche, leurs distances peuvent être calculés.

De nombreux sites web utilisent cette notion de distance afin de proposer par exemple des services.

Calcul de distance

Le langage Java est utilisé pour les méthodes de calcul. Il peut être déporté au niveau de la base de données en tant que requête ou bien en procédure stockée.

Méthode de calcul

Il y a différente façon de faire le calcul selon la précision et les performances que l'on désire.

Formule d'Haversine

Cette méthode est assez précise.

a = sin²(Δφ/2) + cos(φ1).cos(φ2).sin²(Δλ/2)
c = 2.atan2(√a, √(1−a))
d = R.c

public void distanceHaversine() {
        double R = 6371; // km
        double lat1 = 45.3999600000;
        double lat2 = 45.4223330000;
        double lon1 = -73.5065150000;
        double lon2 = -73.4887840000;
        double dLat = Math.PI * (lat2 - lat1) / 180;
        double dLon = Math.PI * (lon2 - lon1)/ 180;
        double lat1Rad = Math.PI * lat1 / 180;
        double lat2Rad = Math.PI * lat2 / 180;
        double a = Math.sin(dLat / 2) * Math.sin(dLat / 2)
                + Math.sin(dLon / 2) * Math.sin(dLon / 2) * Math.cos(lat1Rad) * Math.cos(lat2Rad);
        double c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a));
        double d = R * c;
        System.out.println(d + " Km ");
}

Le résultat obtenu est: 2.8468710743613688 Km 

Loi des cosinus sphériques

d = acos( sin(φ1).sin(φ2) + cos(φ1).cos(φ2).cos(Δλ) ).R
R est le rayon de la terre soit: 6371

Les valeurs utilisées dans cette formule doivent être des radians. Les valeurs de latitude et de longitude du fichier csv doit être convertie. Nous allons reprendre les données présentés ci-haut. 


public void distanceLoisCosinusSphérique(){
        int R = 6371; //KM
        double lat1=45.3999600000;
        double lat2=45.4223330000;
        double lon1=-73.5065150000;
        double lon2=-73.4887840000;
        //conversion des valeurs en radians
        double lat1Rad = Math.radians(lat1);
        double lat2Rad = Math.radians(lat2);
        double lon1Rad = Math.radians(lon1);
        double lon2Rad = Math.radians(lon2);
        //calcule     
        double d = Math.acos(Math.sin(lat1Rad)*Math.sin(lat2Rad) +
                  Math.cos(lat1Rad)*Math.cos(lat2Rad) *
                  Math.cos(lon2Rad-lon1Rad)) * R;
        System.out.println(d + " Km ");
}

Le résultat obtenu est: 2.8468710738500866 Km  


Projection cylindrique équidistante

Cette méthode est performante mais pas précise pour de grande distance. Elle peut donc s'avérer intéressante sur de petite distance.

x = Δλ.cos(φ)
y = Δφ
d = R.√x² + y²


public void distanceProjectionCylEqui(){
    int R = 6371; // km
    double lat1=45.3999600000;
    double lat2=45.4223330000;
    double lon1=-73.5065150000;
    double lon2=-73.4887840000;
    //conversion des valeurs en radians
    double lat1Rad = Math.radians(lat1);
    double lat2Rad = Math.radians(lat2);
    double lon1Rad = Math.radians(lon1);
    double lon2Rad = Math.radians(lon2);
    //calcule      
    double x = (lon2Rad - lon1Rad) * Math.cos((lat1Rad + lat2Rad) / 2);
    double y = (lat2Rad - lat1Rad);
    double d = Math.sqrt(x * x + y * y) * R;
    System.out.println(d + " Km ");

}

Le résultat obtenu est:  2.8468710926857486 Km


Distance max autour d'un point (radius)

Cette fonctionnalité permet de trouver tous les codes postaux à une certaine distance d'un autre.
Par exemple tous les codes postaux à 1 kilomètre du mien, toutes les coiffeurs à 2 kilomètres de ma maison.

Le fichier csv de cet article sera inséré dans une bd H2.

CREATE TABLE POSTAL_CODE_CANADA AS SELECT * FROM CSVREAD('Canada.csv');

Exécuter ces commandes sous H2, car il importe les columns en tant que varchar2 ce qui n'est pas optimale pour effectuer des calculs.


ALTER TABLE POSTAL_CODE alter column latitude double
ALTER TABLE POSTAL_CODE alter column LONGITUDE double


Ne pas oublier de modifier la première ligne de ce fichier afin de spécifier les noms de colonnes.Sinon une structure similaire à celle là peut être créée.

CREATE TABLE POSTAL_CODE_CANADA
  POSTAL_CODE VARCHAR (32) NOT NULL,
  LATITUDE DOUBLE NOT NULL,
  LONGITUDE DOUBLE NOT NULL,
  CITY VARCHAR(32),
  PROVINCE VARCHAR(32)
);

Une requête SQL pour rechercher les codes postaux inférieur à 5km d'un point.

SELECT 
FROM 
  SELECT  ostal_code, 6378 * acos( cos(radians(latitude)) * cos(radians(latitudeP)) * cos   (radians(longitudeP)-  radians(longitude)) + sin(radians(latitude)) * sin(radians(latitudeP))) as   Distance 
FROM postal_code_CA 
WHERE distance <=5

Dans cette requête, remplacer latitudeP par lattitude et longitudeP par la longitude du point désiré.

SELECT *
FROM (
   select POSTAL_CODE, CITY, PROVINCE, 6378 * acos( cos(radians(latitude)) * cos(radians(45.39996)) * cos (radians(-73.506515)-radians(longitude)) + sin(radians(latitude)) * sin(radians(45.39996))) as Distance
FROM POSTAL_CODE_CA
) t
WHERE Distance < 5 ORDER BY DISTANCE

Une telle requête sous H2 prenait jusqu'à 3 secondes. Assurez-vous de bien optimisez votre base de données afin d'obtenir des délais raisonnables.

Fonctionnalité de base de donnée

Certaine base de donnée ont des fonctions pour effectuer des calculs encore plus précis et rapide. Par exemple, PostgreSQL a l'extension PostGIS qui utilise les index de type R-Tree.
Oracle a des capacités similaire ainsi que Microsoft SQL Server 2008 (colonne de type Geography).

Aucun commentaire:

Enregistrer un commentaire