social.py 6.56 KB
Newer Older
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')