|
Recherche d'un processus d'historisation de base de données d'occupation des sols appliqué au référentiel géographique forestier de l'IGN( Télécharger le fichier original )par Romain Louvet Université Paris Diderot - Paris 7 - M1 Géographie et Sciences des territoires 2013 |
VI.A - « identite_1.py »# User : Romain Louvet, stagiaire # 28/08/2013 # IGN, equipe produit forêt environnement, projet OCS & mise a jour du RGFor # DESCRIPTION : # test d'algorithme d'appariement attribuant automatiquement une identite pre-existante a un nouvel objet. # Les criteres d'appariement sont : # - CRITERE 1 : le successeur possede la meme nomenclature que l'anterieur # - CRITERE 2 : le successeur est compris dans l'anterieur, le successeur recupere l'identite de l'anterieur # - CRITERE 3 : le successeur contient l'anterieur et l'aire de l'anterieur est > ou = a 40% de l'aire du successeur # - CRITERE 4 : le successeur intersecte l'anterieur et l'aire du successeur est > ou = a 40% de l'aire du anterieur # Si CRITERE 1 ET (CRITERE 2 ou 3 ou 4) sont vrais, l'entite successeur = l'entite anterieur # Sinon, le successeur est une nouvelle entite. # IMPORTANT : # - ouvrir session de mise a jour avant de charger # - charger dans fenetre python d'ArcGis sous arcmap (mode immediat d'utilisation de python avec ArcGis) # - charger indentite_1, activer mode edition, charger intendite_2 # - si erreur avec 'locks' au moment de sauvegarder, fermer arcmap, redemarrer, effacer feature class cree, relancer OU # ouvrir dossier de la gdb rechercher "lock" et supprimer lock sur "APRES..." et "AVANT..." (solution non teste) import arcpy from datetime import datetime import sys ### I - VARIABLES a definir fc = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/RGFOR65_test" # adresse de la table d'actualites date = """ date '2010-08-01' """# date de la PVA utilisee pour la mise a jour date_fmt = '%Y-%m-%d' gdb = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/" # adresse d'enregistrement AVANT = "AVANT2010" # nom de la table des objets valide jusqu'a la date de PVA APRES = "APRES2010" # nom de la table des objets valide depuis la date de PVA ### II - creation couche AVANT - APRES date de la PVA a partir d'une selection arcpy.MakeFeatureLayer_management (fc, "actulyr") x arcpy.SelectLayerByAttribute_management ("actulyr", "NEW_SELECTION", """ "DATE_V_TO" = """+date) arcpy.CopyFeatures_management ("actulyr",gdb+AVANT) arcpy.SelectLayerByAttribute_management ("actulyr", "NEW_SELECTION", """ "DATE_V_FROM" = """+date) arcpy.CopyFeatures_management ("actulyr",gdb+APRES) identite_2.py identite_3.py evenements_test.py xi VI.B - « identite_2.py » import arcpy from datetime import datetime import sys ### I - VARIABLES a definir fc = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/RGFOR65_test" # adresse de la table d'actualites date = '2010-08-01'# date de la PVA utilisee pour la mise a jour date_fmt = '%Y-%m-%d' gdb = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/" # adresse d'enregistrement AVANT = "AVANT2010" # nom de la table des objets valide jusqu'a la date de PVA APRES = "APRES2010" # nom de la table des objets valide depuis la date de PVA # liste des codes d'entites code_ent = open(" C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/code_RGFOR065.txt","r") code_ent = code_ent.read() code_ent = code_ent.split("\n") ### III - selection par localisation dans AVANT depuis APRES # creation d'une liste contenant les identifiants absolus des lignes listeOID = [] # creation d'un curseur de recherche pour parcourir la table "APRES" lignes = arcpy.SearchCursor(APRES) # boucle recherche ligne par ligne de la table "APRES" et ajout des valeurs a la liste "listeOID" for ligne in lignes: # recupere la valeur du champ OBJECTID, champ d'identifiant absolu value = ligne.getValue("OBJECTID") # ajoute a la liste "listeOID" listeOID.append(value) # fait de meme avec "CLEABS" value = ligne.getValue("CLEABS") # mesure "listeOID" pour determiner le nombre d'iteration lgr = len(listeOID) # trouve la derniere valeur "ENTITE_ID" # servira a appeler nouveau code pour nouvelle entite arcpy.SelectLayerByAttribute_management ("actulyr", "NEW_SELECTION", """ "ENTITE_ID" = (SELECT MAX ("ENTITE_ID") FROM RGFOR65_test) """) resultat = arcpy.SearchCursor("actulyr") for cmax in resultat: cle = cmax.getValue("OBJECTID") xii # boucle for i in range(lgr): # selection par attribut dans la table APRES de la ligne d'index "i" dans "listeOID" OID = str(listeOID[i]) select = """ "OBJECTID" = """+OID arcpy.SelectLayerByAttribute_management (APRES, "NEW_SELECTION", select) # creation d'un curseur de recherche dans la table "APRES" lignes = arcpy.SearchCursor(APRES) # condition si l'objet est une nouvelle entite new = 0 # boucle de recherche ligne par lignee de la table "APRES" des valeurs "NOMEN", "SHAPE_Area" et "DATE_V_CREA" for ligne in lignes: tfvAP = ligne.getValue("NOMEN") aireAP = ligne.getValue("SHAPE_Area") # CRITERES 3 et 4 : calcul du % de l'aire dans "APRES" aireAP100a = aireAP*40/100 aireAP100b = aireAP/1.4 # CRITERE 2 : selection par emplacement dans la couche "AVANT" des polygones # qui contiennent polygones dans la couche "APRES" arcpy.SelectLayerByLocation_management (AVANT, "CONTAINS", APRES) # creation d'un curseur pour la couche "AVANT" lignes2 = arcpy.SearchCursor(AVANT) # condition si selection ne fonctionne pas entIDAV = 0 # pour toutes les lignes de la selection for ligne2 in lignes2: # recupere valeur du champ "NOMEN" et "DATE_CREA" tfvAV = ligne2.getValue("NOMEN") # CRITERE 1 : si valeur "NOMEN" est identique a celle de l'objet dans "APRES" if tfvAP == tfvAV: # recupere valeur du champ "ENTITE_ID", "PARTIE" et "DATE_V_CREA" de "AVANT" et l'attribue a celui de "APRES" entIDAV = ligne2.getValue("ENTITE_ID") entIDAV = str(entIDAV) entiteID = '"'+""+entIDAV+""+'"' partieAV = ligne2.getValue("PARTIE") partieAV = str(partieAV) partie = '"'+""+partieAV+""+'"' crea = ligne2.getValue("DATE_V_CREA") xiii crea = '"'+""+str(crea)+""+'"' arcpy.CalculateField_management (APRES,"ENTITE_IDtest", entiteID, "PYTHON") arcpy.CalculateField_management (APRES,"PARTIE_test", partie, "PYTHON") arcpy.CalculateField_management (APRES,"DATE_V_CREA_test", crea, "PYTHON") # l'objet n'est pas une nouvelle entite new = 1 # si CONTAINS ne donne rien : if entIDAV == 0: # CRITERE 3 : selection par emplacement dans la couche "AVANT" des polygones # qui sont contenus par des polygones dans la couche "APRES" arcpy.SelectLayerByLocation_management (AVANT, "WITHIN", APRES) # creation d'un curseur pour la couche "AVANT" lignes2 = arcpy.SearchCursor(AVANT) # condition si selection ne fonctionne pas entIDAV = 0 # listes pour la selection liAIREAV = [] liENTITEAV = [] liPARTIE = [] liCREA = [] # pour toutes les lignes dans la selection for ligne2 in lignes2: # recupere valeur du champ "NOMEN" tfvAV = ligne2.getValue("NOMEN") # CRITERE 1 : si valeur "NOMEN" de "AVANT" est identique a celle de l'objet # dans "APRES" ajoute valeur de son champ "SHAPE_Area" et "ENTITE_ID" a une liste if tfvAP == tfvAV: aireAV = ligne2.getValue("SHAPE_Area") entIDAV = ligne2.getValue("ENTITE_ID") partieAV = ligne2.getValue("PARTIE") crea = ligne2.getValue("DATE_V_CREA") liAIREAV.append(aireAV) liENTITEAV.append(entIDAV) liPARTIE.append(partieAV) liCREA.append(crea) # trouve l'objet le plus grand et son "ENTITE_ID" try: airemax = 0 airemax = max(liAIREAV) xiv ind = liAIREAV.index(airemax) except (SyntaxError,ValueError,IndexError): pass # si aire "AVANT" remplie condition par rapport a aire "APRES" if airemax > aireAP100a: # recupere valeur du champ "ENTITE_ID", "PARTIE" et "DATE_V_CREA" de "AVANT" et l'attribue a celui de "APRES" entIDAV = liENTITEAV[ind] partieAV = liPARTIE[ind] crea = liCREA[ind] entIDAV = str(entIDAV) entiteID = '"'+""+entIDAV+""+'"' partieAV = str(partieAV) partie = '"'+""+partieAV+""+'"' crea = '"'+""+str(crea)+""+'"' arcpy.CalculateField_management (APRES,"ENTITE_IDtest", entiteID, "PYTHON") arcpy.CalculateField_management (APRES,"PARTIE_test", partie, "PYTHON") arcpy.CalculateField_management (APRES,"DATE_V_CREA_test", crea, "PYTHON") # l'objet n'est pas une nouvelle entite new = 1 # si WITHIN ne donne rien non plus if entIDAV == 0: # CRITERE 4 : selection par emplacement dans la couche "AVANT" des polygones # qui possedent une intersection avec les polygones dans la couche "APRES" arcpy.SelectLayerByLocation_management (AVANT, "INTERSECT", APRES) # creation d'un curseur pour la couche "AVANT" lignes2 = arcpy.SearchCursor(AVANT) # listes pour la selection liAIREAV = [] liENTITEAV = [] liPARTIE = [] liCREA = [] for ligne2 in lignes2: # recupere valeur du champ "NOMEN" tfvAV = ligne2.getValue("NOMEN") # CRITERE 1 : si valeur "TFV" de "AVANT" est identique a celle de l'objet # dans "APRES" ajoute valeur de son champ "SHAPE_Area" et "ENTITE_ID" a une liste if tfvAP == tfvAV: aireAV = ligne2.getValue("SHAPE_Area") entIDAV = ligne2.getValue("ENTITE_ID") xv partieAV = ligne2.getValue("PARTIE") crea = ligne2.getValue("DATE_V_CREA") liAIREAV.append(aireAV) liENTITEAV.append(entIDAV) liPARTIE.append(partieAV) liCREA.append(crea) # trouve l'objet le plus grand et son "ENTITE_ID" try: airemax = 0 airemax = max(liAIREAV) ind = liAIREAV.index(airemax) except (SyntaxError,ValueError,IndexError): pass # si aire "AVANT" remplie condition par rapport a aire "APRES" if airemax >= aireAP100b: # recupere valeur du champ "ENTITE_ID", "PARTIE" et "DATE_V_CREA" de "AVANT" et l'attribue a celui de "APRES" entIDAV = liENTITEAV[ind] partieAV = liPARTIE[ind] crea = liCREA[ind] entIDAV = str(entIDAV) entiteID = '"'+""+entIDAV+""+'"' partieAV = str(partieAV) partie = '"'+""+partieAV+""+'"' crea = '"'+""+str(crea)+""+'"' arcpy.CalculateField_management (APRES,"ENTITE_IDtest", entiteID, "PYTHON") arcpy.CalculateField_management (APRES,"PARTIE_test", partie, "PYTHON") arcpy.CalculateField_management (APRES,"DATE_V_CREA_test", crea, "PYTHON") # l'objet n'est pas une nouvelle entite new = 1 # si l'objet est une nouvelle entite if new == 0: # cherche un nouveau code ENTITE_ID dans la liste entiteID = '"'+""+str(code_ent[cle])+""+'"' # attribue premier element de la liste de code pour PARTIE partie = '"'+""+"aaa"+""+'"' # attribue nouvelle date de creation crea = datetime.strptime(date,date_fmt) crea = '"'+""+str(crea)+""+'"' xvi arcpy.CalculateField_management (APRES,"ENTITE_IDtest", entiteID, "PYTHON") arcpy.CalculateField_management (APRES,"PARTIE_test", partie, "PYTHON") arcpy.CalculateField_management (APRES,"DATE_V_CREA_test", crea, "PYTHON") # prochain code cle = cle + 1 xvii VI.C - « identite_3.py » import arcpy from datetime import datetime import sys ### I - VARIABLES a definir fc = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/RGFOR65_test" # adresse de la table d'actualites date = '2010-08-01'# date de la PVA utilisee pour la mise a jour date_fmt = '%Y-%m-%d' gdb = " C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/test_rgfor65.gdb/test65/" # adresse d'enregistrement AVANT = "AVANT2010" # nom de la table des objets valide jusqu'a la date de PVA APRES = "APRES2010" # nom de la table des objets valide depuis la date de PVA ### IV- Partie ### Si plusieurs objets ont le meme code ENTITE_ID, ils doivent avoir differents codes PARTIE ### Variables # liste des codes PARTIE code_partie = open(" C:/Users/Romain/Desktop/IGN/test_histoRGFOR_OCS/code_PARTIE065.txt","r") code_partie = code_partie.read() code_partie = code_partie.split("\n") # Variables pour IVA liOID = [] # liste OBJECTID liENT = [] # liste ENTITE liPAR = [] # liste PARTIE depart liENPAR = [] # liste id d'entite complet liARE = [] # liste SHAPE_area # Variable pour IVB liTTL = [] # liste totale index liPAR2 = [] # liste PARTIE arrivee ### IVA - Listes des valeurs des champs # creation d'un curseur de recherche pour parcourir la table "APRES" lignes = arcpy.SearchCursor(APRES) # boucle recherche ligne par ligne de la table "APRES" for ligne in lignes: # recupere la valeur du champ "OBJECTID", "ENTITE_IDtest", "PARTIE_test" et "SHAPE_Area" vOID = ligne.getValue("OBJECTID") vENT = ligne.getValue("ENTITE_IDtest") vPAR = ligne.getValue("PARTIE_test") vARE = ligne.getValue("SHAPE_Area") # ajoute aux listes xviii liOID.append(vOID) liENT.append(vENT) liPAR.append(vPAR) liENPAR.append(vENT+vPAR) liARE.append(vARE) liPAR2 = liPAR lgr = len(liOID) # longueur de la table ### IVB - Boucle de detection de nouvelle partie for i in range(lgr): # creation et remise a zero des variables calculees par la boucle liIND = [] # liste index des lignes avec meme ENTITE_ID + PARTIE liARE2 = [] # valeur aire de ces lignes liPAR1 = [] # liste PARTIE en cours # si cette ligne n'a pas deja ete traitee # essaye de trouver i dans liste totale index try: test = liTTL.index(i) # si ne fonctionne pas, lance la recherche except(ValueError): # cherche si code ENTITE_ID plus PARTIE # existe plusieurs fois enpar = liENPAR[i] # pour chacun des codes ENTITE_ID + PARTIE for indx, ID in enumerate(liENPAR): # si le code est egal au code recherche if ID == enpar: # son index est ajoute a la liste liIND.append(indx) # si cette liste est > 1, il existe plusieurs codes identique lgi = len(liIND) if lgi > 1: # avec la liste d'index je trouve la valeur max du code PARTIE et de l'airepour cette ENTITE_ID for j in liIND: liARE2.append(liARE[j]) liPAR1.append(liPAR[j]) mxaire = max(liARE2) mxliPAR1 = max(liPAR1) # trouve valeur PARTIE la plus grande dans liste num = code_partie.index(mxliPAR1) xix # pour toutes les lignes, sauf celle avec la plus grande aire # le champ PARTIE possede un nouveau code for k in liIND: aire = liARE[k] if aire < mxaire: # donne index du code suivant dans la liste de code PARTIE num = num+1 # remplace ancien par nouveau code liPAR2[k]= code_partie[num] # ajoute les lignes deja faites liTTL.extend(liIND) ### IVC - boucle calcul du champ PARTIE for i in range(lgr): arr = liPAR2[i] # selection de la ligne a modifier OID = str(liOID[i]) select = """ "OBJECTID" = """+OID arcpy.SelectLayerByAttribute_management (APRES, "NEW_SELECTION", select) lignes = arcpy.SearchCursor(APRES) # modification du champ "PARTIE" avec le nouveau code arr = '"'+""+str(arr)+""+'"' arcpy.CalculateField_management (APRES,"PARTIE_test", arr, "PYTHON") xx |
|
