Code TD3 géocodage

# -*- coding: utf-8 -*-

__author__ = 'donald'

import sys, os
sys.path.append(r'C:\Program Files (x86)\QGIS Lyon\apps\Python27\Lib\site-packages')
os.environ['PATH'] = r'C:\Program Files (x86)\QGIS Lyon\bin'
os.environ['GDAL_DATA']=r'C:\Program Files (x86)\QGIS Lyon\share\gdal'


from osgeo import ogr,osr

# address_no=float(raw_input("Saisir le numero de l'adresse:"))
# address_street=raw_input("Saisir la rue de l'adresse:")
# address_street1=address_street.replace(" ","")
# address_city=raw_input("Saisir la ville de l'adresse:")
# address_city1=address_city.replace(" ","")

StreetFile = r'C:\Users\donald\Documents\Université Sherbrooke\Travail session Python\RRN_QC_8_0_SEGMROUT_MTL.shp'
# transformation de système de projection
source = osr.SpatialReference()
source.ImportFromEPSG(4326)#lat long
target = osr.SpatialReference()
target.ImportFromEPSG(2154)# en mètres
transform = osr.CoordinateTransformation(source, target)

dataSource1 = ogr.Open(StreetFile)
daLayer1 = dataSource1.GetLayer()

# layerDefinition = daLayer1.GetLayerDefn()
# for i in range(layerDefinition.GetFieldCount()):
#     print layerDefinition.GetFieldDefn(i).GetName()
# print ""


point1=ogr.Geometry(ogr.wkbPoint)
latitude=45.513522
longitude=-73.605681
point1.AddPoint (longitude,latitude)#(-73.602473,45.507056)#(-73.561060,45.512703)#(-73.554,45.508)45.512703,
point1.Transform(transform)
print point1
pt_txt=str(point1)
pt_txt1=pt_txt.replace("POINT (","")
pt_txt2=pt_txt1.replace(","," ")
pt_txt3=pt_txt2.replace(")","")
list_pt=pt_txt3.split()
print list_pt

bufferDistance = 200
buffer_geom = point1.Buffer(bufferDistance)


# pour trouver les rues situées dans le buffer
min_distance=bufferDistance
street_ID=0
for feature in daLayer1:
    pt_geom=feature.GetGeometryRef()
    pt_geom.Transform(transform)
    if pt_geom.Within(buffer_geom):
        print "ID:",feature.IDN
        print feature.NOMRUE_C_G
        print pt_geom.Length()
        print pt_geom
        print "Distance:",pt_geom.Distance(point1)
        print ""
        if pt_geom.Distance(point1)<min_distance:
            min_distance=pt_geom.Distance(point1)
            street_ID=feature.IDN
print "Smallest distance:",min_distance,"ID:",street_ID

# rouvrir de nouveau le layer pour faire la boucle et trouver la rue la plus proche
dataSource1 = ogr.Open(StreetFile)
daLayer1 = dataSource1.GetLayer()

# pour trouver le segment le plus proche
for item in daLayer1:
    if item.IDN==street_ID:
        print "Closest street:",item.NOMRUE_C_G
        geom=item.GetGeometryRef()
        geom.Transform(transform)
        line_txt=str(geom)
        line_txt1=line_txt.replace("LINESTRING (","")
        line_txt2=line_txt1.replace(","," ")
        line_txt3=line_txt2.replace(")","")
        list=line_txt3.split()
        tot_length=geom.Length()
        print list
        print "Length:",tot_length

        #trouver le segment approprié par rapport au point demandé
        spread=float(feature.NUMD_G-feature.NUMP_G)
        #if spread!=0:
        # target_length=((address_no-feature.NUMP_G)/spread)*tot_length
        # print target_length
        j=0
        cum_long=0
        while j<(len(list)-2):
                #longueur du segment
                # long=((X2-X1)**2+(Y2-Y1)**2)**0.5
                long=((float(list[j+2])-float(list[j]))**2+(float(list[j+3])-float(list[j+1]))**2)**0.5
                print "Length of segment:",long
                #distance du point au segment
                #   var=X2-X1
                var=float(list[j+2])-float(list[j])
                # print "X2-X1:",var
                # print "Y2-Y1/var:",(float(list[j+3])-float(list[j+1]))/var
                # print "X3:",float(list_pt[1])
                # print "X1:",float(list[j])
                #dist=(Y3-Y1-[((Y2-Y1)/var)*(X3-X1))*var/long
                dist=(float(list_pt[1])-float(list[j+1])-(((float(list[j+3])-float(list[j+1]))/var)*(float(list_pt[0])-float(list[j]))))*var/long
                if dist<0:
                    dist=dist*(-1)
                print "Distance to segment:",dist
                print ""
                j+=2