social.py 6.56 KB
 Seneque Colin committed Feb 28, 2019 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 ``````#!/usr/bin/python3 # -*- coding: utf-8 -*- from shapely.geometry import mapping, shape from shapely.ops import cascaded_union import fiona from fiona import collection import rtree import math import geopandas as gpd # Transition matrix : cost of transformation between patch X and patch Y # 0: cereales, 1: cultures industrielles, 2: Fourrages, 3: fruits_legumes, 4: oleagineux, 5: prairies, 6: proteagineux matrice_transition = {} for i in range(7): matrice_transition[i]={} # Céréales matrice_transition[0][0] = 0 matrice_transition[0][1] = 0.5 matrice_transition[0][2] = 0.5 matrice_transition[0][3] = 1 matrice_transition[0][4] = 0.5 matrice_transition[0][5] = 0 matrice_transition[0][6] = 0.5 # Cultures industrielles matrice_transition[1][0] = 0.5 matrice_transition[1][1] = 0 matrice_transition[1][2] = 0.5 matrice_transition[1][3] = 1 matrice_transition[1][4] = 0.5 matrice_transition[1][5] = 0 matrice_transition[1][6] = 0.5 # Fourrages matrice_transition[2][0] = 0.5 matrice_transition[2][1] = 0.5 matrice_transition[2][2] = 0 matrice_transition[2][3] = 1 matrice_transition[2][4] = 0.5 matrice_transition[2][5] = 0 matrice_transition[2][6] = 0.5 # Fruits et légumes matrice_transition[3][0] = 0.5 matrice_transition[3][1] = 0.5 matrice_transition[3][2] = 0.5 matrice_transition[3][3] = 0 matrice_transition[3][4] = 0.5 matrice_transition[3][5] = 0 matrice_transition[3][6] = 0.5 # Oléagineux matrice_transition[4][0] = 0.5 matrice_transition[4][1] = 0.5 matrice_transition[4][2] = 0.5 matrice_transition[4][3] = 1 matrice_transition[4][4] = 0 matrice_transition[4][5] = 0 matrice_transition[4][6] = 0.5 # Prairies matrice_transition[5][0] = 0.5 matrice_transition[5][1] = 0.5 matrice_transition[5][2] = 0.5 matrice_transition[5][3] = 1 matrice_transition[5][4] = 0.5 matrice_transition[5][5] = 0 matrice_transition[5][6] = 0.5 # Protéagineux matrice_transition[6][0] = 0.5 matrice_transition[6][1] = 0.5 matrice_transition[6][2] = 0.5 matrice_transition[6][3] = 1 matrice_transition[6][4] = 0.5 matrice_transition[6][5] = 0 matrice_transition[6][6] = 0 # Transform the culture's name in his corresponding code def return_culture_code(culture): if culture == "Cereales" : #print culture return 0 if culture == "Cultures industrielles" : return 1 if culture == "Fourrages" : return 2 if culture == "Fruits et legumes" : #print culture return 3 if culture == "Oleagineux" : #print culture return 4 if culture == "Prairies" : return 5 if culture == "Proteagineux" : return 6 def social(parcelle_initiale_expl, scenario, bdv): ''' This function is used to calculate the social indice. The idea is to propose a transition matrix for each type of culture and use it to visualize the conversion cost of the total conversion between the initial scenario and the old one. Each transition as a value between 0 and 1. These values were defined during reunions and with the review of our experts. :param parcelle_initiale_expl: ShapeFile of the PAT with the ID of exploitant as an attribute for the patches :param scenario: The scenario to analyse as a list :param bdv: ShapeFile of the living area (bassin de vie) ''' with fiona.open(bdv, 'r') as layer1: with fiona.open(parcelle_initiale_expl, 'r' ) as layer_parcelle_initiale : list_parcelle_initiale = [] for row in layer_parcelle_initiale : list_parcelle_initiale.append(row) list_parcelle_initiale = sorted(list_parcelle_initiale, key=lambda k: k['properties']['ID_PARCEL']) list_culture = [0,1,2,3,4,5,6] # Contains all kind of patches # List_expl is the list who will contains each exploitant of the PAT and his patches list_expl = {} for i in range(len(list_parcelle_initiale)) : id_expl = list_parcelle_initiale[i]["properties"]["id_expl"] if id_expl not in list_expl: list_expl[id_expl] = [] list_expl[id_expl].append(i) # culture_bdv will contains all the details about all kinds of cultures in each of the bdv in the initial PAT, area for each culture + total area of the bassin. # It will be easy to get % of each culture in each bdv then. culture_bdv = {} index = rtree.index.Index() for feat1 in layer1: # For each bdv fid = int(feat1['id']) # In the dictionnary 'resilience', create a field 'surf_totale' who will contains the sum of area of each kind of culture in this square culture_bdv[feat1["properties"]['BV2012BV20']] = {} culture_bdv[feat1["properties"]['BV2012BV20']]['surf_totale'] = 0.0 # This loop add the area of each kind of culture in a specific field for each square of the grid for c in list_culture : culture_bdv[feat1["properties"]['BV2012BV20']][c] = 0.0 for feat2 in list_parcelle_initiale: # for each element of the scenario 'list_parcelle_initiale if feat2["properties"]["rpg_lfgc_s"] != None : culture_bdv[feat2["properties"]['Bdv']][return_culture_code(feat2['properties']['rpg_lfgc_s'])] += float(feat2['properties']['SURF_PARC']) culture_bdv[feat2["properties"]['Bdv']]['surf_totale'] += float(feat2['properties']['SURF_PARC']) i = 0 cpt = 0 somme_indice = 0 value_indice = 0 for patches in scenario : value_indice = 0 if list_parcelle_initiale[i]['properties']['rpg_lfgc_s'] != None : code_culture_initial = return_culture_code(list_parcelle_initiale[i]['properties']['rpg_lfgc_s']) code_culture_scenario = return_culture_code(scenario[i]['properties']['rpg_lfgc_s']) if code_culture_initial != code_culture_scenario : cpt += 1 value_indice = matrice_transition[code_culture_initial][code_culture_scenario] id_expl_patch = list_parcelle_initiale[i]["properties"]["id_expl"] for fid in list_expl[id_expl_patch] : if list_parcelle_initiale[fid]["properties"]["rpg_lfgc_s"] != None : if return_culture_code(list_parcelle_initiale[fid]["properties"]["rpg_lfgc_s"]) == code_culture_scenario : value_indice = 0 break bdv_patch = scenario[i]["properties"]["Bdv"] # In the case were there is already a lot of culture of a kind in a living area, the cost from any culture to this one will be 0. if culture_bdv[str(bdv_patch)][code_culture_scenario] / culture_bdv[str(bdv_patch)]["surf_totale"] > 0.1 : value_indice = 0 somme_indice += value_indice i += 1 if cpt != 0 : return somme_indice/cpt else : return 0 #social('./data/Parcelle_PAT_exploitant.shp', '../../Etape1_shp_propre/MCMC/testSimulationSurf_MCMC.shp', './data/Bassin_de_vie_Pat.shp')``````