-
Notifications
You must be signed in to change notification settings - Fork 2
/
03_3_graphs_do_maps.R
1214 lines (980 loc) · 64.2 KB
/
03_3_graphs_do_maps.R
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
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
# Ken <- readRDS('//dapadfs/workspace_cluster_8/climateriskprofiles/results/all_countrys_maps/Kenya_P/gadm36_KEN_2_sp.rds')
# Maps for
# H. Achicanoy & A. Esquivel
# Alliance Bioversity-CIAT
# July - 2020.
rm(list = ls()); gc(reset = TRUE)
# =--------------------
# Packages
options(warn = -1, scipen = 999)
suppressMessages(library(pacman))
suppressMessages(pacman::p_load(tidyr, dplyr, tibble, ggplot2, raster, ncdf4, sf, lubridate, glue, cowsay, fst, ggspatial, vroom, sp, compiler))
# =--------------------
# =----------------------------------
# Identificacion de pixel para ETH
# =----------------------------------
country <- 'Ivory_Coast'
count_i <- c('Lagunes' ,
'Comoe' )
iso3c <- 'CIV'
Big <- 'B'
adm_lvl <- 1
for(i in 1:length(count_i)){
county <- count_i[i]
# -----------------------------------
C_shp <- county
co <- tolower(country)
country1 <- co
# =--------------------
# Ruta Principal para guardados:
root <- '//dapadfs/workspace_cluster_8/climateriskprofiles/'
# Load county shapefile
if(country == 'India'){
# India
country1 <- raster::shapefile(glue::glue('//dapadfs/workspace_cluster_8/climateriskprofiles/results/India/states/Admin2.shp'))
shp <- raster::shapefile(glue::glue('//dapadfs/workspace_cluster_8/climateriskprofiles/results/India/states/Admin2.shp'))
shp <- shp[shp@data$ST_NM %in% C_shp,]
plot(shp)
shp@data$ISO <- iso3c
}else{
country1 <- readRDS(glue::glue('{root}data/shps/shps_from_R/{country}/gadm36_{iso3c}_{adm_lvl}_sp.rds'))
shp <- readRDS(glue::glue('{root}data/shps/shps_from_R/{country}/gadm36_{iso3c}_{adm_lvl}_sp.rds'))
shp@data$NAME_1 <- iconv(shp@data$NAME_1,from="UTF-8",to="ASCII//TRANSLIT")
shp@data$NAME_1 <- case_when(shp@data$NAME_1 == 'Trans Nzoia' ~ 'Trans-Nzoia',
shp@data$NAME_1 == "Murang'a" ~ 'Muranga',
TRUE ~ shp@data$NAME_1)
if(iso3c == 'ETH'){country1$NAME_1 = country1$NAME_2; shp$NAME_1 = shp$NAME_2}else{print('ok')}
shp <- shp[shp@data$NAME_1 %in% C_shp,]
plot(shp)
shp@data$ISO <- iso3c
}
# Load id coords
crd <- vroom('//dapadfs/workspace_cluster_8/climateriskprofiles/data/id_all_country.csv')
crd <- crd %>%
dplyr::filter(Country == country)
pnt <- crd %>% dplyr::select(x,y) %>% sp::SpatialPoints(coords = .)
crs(pnt) <- crs(shp)
# Filter coordinates that are present in the county
pnt <- sp::over(pnt, shp) %>% data.frame %>% dplyr::select(ISO) %>% complete.cases() %>% which()
crd <- crd[pnt,]
crd <<- crd
# =------------------------------------------
glwd1 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_1.shp' )
crs(glwd1) <- crs(shp)
glwd1 <- rgeos::gSimplify(glwd1, tol = 0.05, topologyPreserve = TRUE) %>%
sf::st_as_sf()
glwd1 <<- glwd1
glwd2 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_2.shp' )
crs(glwd2) <- crs(shp)
glwd2 <- rgeos::gSimplify(glwd2, tol = 0.05, topologyPreserve = TRUE) %>%
sf::st_as_sf()
glwd2 <<- glwd2
# =------------------------------------------
##### =-----------------------------------------------------------
# ### =-----------------------------------------------------------
# =------------ Funciones...
do_climate_maps <- function(data_split){
ISO3 <- unique(data_split$ISO3); county <- unique(data_split$county)
semester <- unique(data_split$semester) ; Country <- unique(data_split$Country)
# Conditions...
if(Big =='N'){
median_data <- data_split %>%
dplyr::group_by(id, county, Country, x, y, ISO3, semester, time) %>%
dplyr::select(-year) %>%
dplyr::summarise_all(mean) %>% dplyr::ungroup() %>%
dplyr::mutate(NT35 = round(x = NT35, digits = 0),
ndws = round(ndws,0))
}else if(Big == 'B'){
# =---------------------------------------------------------
median_data <- data_split
} else{print('Change big argument... >.<')}
median_data <- median_data %>% dplyr::filter(time == 'future') %>%
group_by( id, ISO3, county, Country, x, y, time, semester) %>%
summarise_all(mean) %>% mutate_at(.vars = vars(CDD, NT35), .funs = function(x){round(x, 0)}) %>%
ungroup() %>%
dplyr::rename('CDD_f' = 'CDD', 'P5D_f' = 'P5D' , 'P95_f'= 'P95' , 'NT35_f' = 'NT35',
'ndws_f' = 'ndws' ) %>%
dplyr::select(-time) %>%
dplyr::inner_join(dplyr::filter(median_data , time == 'past') %>% dplyr::select(-time), . ) %>%
dplyr::mutate(CDD_c = CDD_f - CDD, P5D_c = P5D_f - P5D, P95_c = P95_f - P95, NT35_c = NT35_f -NT35, ndws_c = ndws_f - ndws)
limits <- median_data %>%
dplyr::select(id, contains('_f')) %>%
setNames(c('id', 'CDD', 'P5D', 'P95', 'NT35', 'ndws')) %>%
mutate(time = 'f') %>%
bind_rows(dplyr::select(median_data, id, CDD:ndws) %>% mutate(time = 'p') , .) %>%
dplyr::select(-id) %>%
dplyr::summarise_all(.funs = c('min', 'max'))
# Aqui se hace solo la figura base...
shp_sf <- shp %>% sf::st_as_sf()
country <- country1 %>% sf::st_as_sf()
xlims <- sf::st_bbox(shp_sf)[c(1, 3)]
ylims <- sf::st_bbox(shp_sf)[c(2, 4)]
limx <- sf::st_bbox(country)[c(1, 3)]
limy <- sf::st_bbox(country)[c(2, 4)]
map_world <- raster::shapefile(glue::glue('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/all_country/all_countries.shp')) %>%
sf::st_as_sf()
# =-------------------------------------------------------------------------------
pos <- which(map_world$ISO3 == iso3c)
map_world <- map_world %>% filter(CONTINENT == map_world[pos, ]$CONTINENT)
pos <- which(map_world$ISO3 == iso3c)
spare_mtx <- st_intersects(map_world, map_world)
all_int <- spare_mtx[pos][[1]]
map_world <- map_world[all_int, ]
proof <- as_tibble(st_centroid(map_world) %>% st_coordinates()) %>%
mutate(iso = map_world$ISO3, name = iconv(map_world$NAME,from="UTF-8",to="ASCII//TRANSLIT")) %>%
mutate(Initals = substr(name, start = 1, stop = 3))
if(ISO3 == 'IND'){
af <- as_tibble(st_centroid(shp_sf) %>% st_coordinates()) %>%
mutate( name = shp_sf$ST_NM) %>%
mutate(Initals = substr(name, start = 1, stop = 3))
}else{
af <- as_tibble(st_centroid(shp_sf) %>% st_coordinates()) %>%
mutate( name = shp_sf$NAME_1) %>%
mutate(Initals = substr(name, start = 1, stop = 3))
}
if(iso3c == 'IND'){
country <- country %>% mutate(NAME_0 = iso3c) %>% group_by(NAME_0) %>% summarise()
}else{
country <- country %>% group_by(NAME_0) %>% summarise()
}
b <- ggplot() + geom_sf(data = map_world, fill = NA, color = gray(.8)) +
geom_sf(data = country, fill = 'lightgray', color = gray(.1), alpha = 0.2) +
geom_sf(data = shp_sf, fill = 'gray' , color = gray(.1)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
coord_sf(xlim = limx, ylim = limy) +
ggrepel::geom_label_repel(data=proof, aes(x=X, y=Y, label=name),
arrow = arrow(length = unit(0.03, "npc"), type = "closed", ends = "first"),
force = 10,
size = 6) +
geom_text(data = af, aes(X, Y, label = Initals, shape = Initals), colour ='black', show.legend = TRUE, size=6) +
ggspatial::annotation_scale(location = "tl", width_hint = 0.5, pad_y = unit(0.3, "in")) +
ggspatial::annotation_north_arrow(location = "tl", which_north = "true",
pad_x = unit(0.1, "in"), pad_y = unit(0.5, "in"), # 0.2 # 0.3
style = north_arrow_fancy_orienteering) +
scale_shape_manual(values = 1:nrow(af),
name=NULL,
labels= paste0(af$name, '(' , af$Initals , ')'))+
labs(x = NULL, y = NULL) +
theme_bw() +
theme(legend.position = 'bottom', text = element_text(size=18),
legend.text = element_text(size=18),
legend.title=element_text(size=18)) + guides(shape = guide_legend(ncol = 3))
ggsave(glue::glue('{path}{Country}/graphs/{tolower(county)}/maps/{county}.png') , width = 8, height = 5.5, dpi = 300)
#===---------------------------------------------------------
#=------------------------------------------------------------
my_limits <- c(round(limits$CDD_min,2)-0.1, round(limits$CDD_max,2)+0.1)
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
index_a <- 'CDD'
a <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = CDD)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_a}\n(days) '), title = 'Historic', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20,
label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_a}_past_S{semester}.png') , width = 10, height = 10, dpi = 300)
# =- Lo mismo para futuro...
a1 <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = CDD_f)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_a}\n(days) '), title = 'Future',x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20,
label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_a}_future_S{semester}.png') , width = 10, height = 10, dpi = 300)
my_limits <- c(round(min(median_data$CDD_c),2)-0.1, round(max(median_data$CDD_c), 2)+0.1)
my_breaks <-round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
a_d <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = CDD_c)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_a}\n(days) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(low = '#000099', mid = 'white', high = '#A50026',
limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_a}_S{semester}.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_a}_{semester}.png'), width=28,height=10,units="in", res = 300) # width = 1580, height = 720,
print(gridExtra::grid.arrange(a, a1, a_d, ncol=3))
dev.off()
#===---------------------------------------------------------
#=------------------------------------------------------------
# Siguiente indice...
index_c <- 'P5D'
my_limits <- c(round(limits$P5D_min,2)-0.1, round(limits$P5D_max,2)+0.1)
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 1)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
c <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P5D)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_c}\n(mm) '), title = 'Historic',x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_c}_past_S{semester}.png') , width = 10, height = 10, dpi = 300)
# =- Futuro.
c1 <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P5D_f)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_c}\n(mm) '), title = 'Future', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_c}_future_S{semester}.png') , width = 10, height = 10, dpi = 300)
my_limits <- c(round(min(median_data$P5D_c),2)-0.1, round(max(median_data$P5D_c), 2)+0.1)
my_breaks <-round(seq(my_limits[1], my_limits[2], length.out= 3), 1)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
c_d <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P5D_c)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_c}\n(mm) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(low = '#A50026', mid = 'white', high = '#000099',
limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20,
label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_c}_S{semester}.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_c}_{semester}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(c, c1, c_d, ncol=3))
dev.off()
#===---------------------------------------------------------
index_d <- 'P95'
my_limits <- c(round(limits$P95_min,2)-0.1, round(limits$P95_max,2)+0.1)
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 1)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
d <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P95)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_d}\n(mm) '), title = 'Historic', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_d}_past_S{semester}.png') , width = 10, height = 10, dpi = 300)
# =----
d1 <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P95_f)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_d}\n(mm) '), title = 'Future', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c( limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_d}_future_S{semester}.png') , width = 10, height = 10, dpi = 300)
my_limits <- c(round(min(median_data$P95_c),2)-0.1, round(max(median_data$P95_c), 2)+0.1)
my_breaks <-round(seq(my_limits[1], my_limits[2], length.out= 3), 1)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
d_d <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = P95_c)) + geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_d}\n(mm) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(limits = my_limits,
breaks = my_breaks,
low = '#A50026', mid = 'white', high = '#000099',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_d}_S{semester}.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_d}_{semester}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(d, d1, d_d, ncol=3))
dev.off()
#===---------------------------------------------------------
#=------------------------------------------------------------
index_e <- 'NT35'
my_limits <- c(round(limits$NT35_min,2)-0.1, round(limits$NT35_max,2)+0.1)
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
e <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = NT35)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_e}\n(days) '), title = 'Historic', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_e}_past_S{semester}.png') , width = 10, height = 10)
# =------------
e1 <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = NT35_f)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_e}\n(days) '), title = 'Future', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_e}_future_S{semester}.png') , width = 10, height = 10)
my_limits <- c(round(min(median_data$NT35_c),2)-0.1, round(max(median_data$NT35_c), 2)+0.1)
my_breaks <-round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
e_d <- ggplot() + geom_tile(data = median_data, aes(x = x, y = y, fill = NT35_c)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_e}\n(days) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(limits = my_limits,
breaks = my_breaks,
low = '#000099', mid = 'white', high = '#A50026',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_e}_S{semester}.png') , width = 10, height = 10)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_e}_{semester}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(e, e1, e_d, ncol=3))
dev.off()
#===---------------------------------------------------------
#=------------------------------------------------------------
index_f <- 'ndws'
my_limits <- c(round(limits$ndws_min,2)-0.1, round(limits$ndws_max,2)+0.1)
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
f <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = ndws)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_f}\n(days) '), title = 'Historic', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20,
label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_f}_past_S{semester}.png') , width = 10, height = 10)
# =------------
f1 <- ggplot() +
geom_tile(data = median_data, aes(x = x, y = y, fill = ndws_f)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_f}\n(days) '), title = 'Future', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/{index_e}_future_S{semester}.png') , width = 10, height = 10)
my_limits <- c(round(min(median_data$ndws_c),2)-0.1, round(max(median_data$ndws_c), 2)+0.1)
my_breaks <-round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
f_d <- ggplot() + geom_tile(data = median_data, aes(x = x, y = y, fill = ndws_c )) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = country, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('{index_f}\n(days) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(limits = my_limits,
breaks = my_breaks,
low = '#000099', mid = 'white', high = '#A50026',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_f}_S{semester}.png') , width = 10, height = 10)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/Dif_{index_f}_{semester}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(f, f1, f_d, ncol=3))
dev.off()
return(median_data)}
# =------------
# Funcion for run basic maps.
do_srad_maps <- function(data_split){
ISO3 <- unique(data_split$ISO3); county <- unique(data_split$county)
Country <- unique(data_split$Country)
Big <- unique(data_split$Big)
# Aqui se hace solo la figura base...
shp_sf <- shp %>% sf::st_as_sf()
pais <- country1 %>% sf::st_as_sf()
xlims <- sf::st_bbox(shp_sf)[c(1, 3)]
ylims <- sf::st_bbox(shp_sf)[c(2, 4)]
# # =--------------------------------------------
#
# glwd1 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_1.shp' )
# crs(glwd1) <- crs(shp)
# glwd1 <- rgeos::gSimplify(glwd1, tol = 0.05, topologyPreserve = TRUE) %>%
# sf::st_as_sf()
#
# glwd2 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_2.shp' )
# crs(glwd2) <- crs(shp)
# glwd2 <- rgeos::gSimplify(glwd2, tol = 0.05, topologyPreserve = TRUE) %>%
# sf::st_as_sf()
#
# # =--------------------------------------------
if(Big =='N'){
# gSeason...
gSeason_i <- data_split %>% dplyr::select( time, x, y,year, gSeason) %>%
group_by(time, x, y, year) %>% summarise(gSeason = max(gSeason)) %>%
ungroup() %>% group_by(time, x, y) %>%
summarise(gSeason = round(mean(gSeason, na.rm = TRUE), 0)) %>% ungroup()
# =----------------------------------------------------------------------
# SLGP -- LGP
two_index <- data_split %>% dplyr::select( time, x, y, gSeason, SLGP, LGP) %>%
group_by(time, x, y, gSeason) %>%
summarise_all(.f = function(x){round(mean(x, na.rm = TRUE), 0)}) %>%
# ungroup(year) %>%
arrange(gSeason) %>%
tidyr::nest(data = c('x', 'y','SLGP', 'LGP')) %>%
drop_na() %>% filter(gSeason < 3) %>% tidyr::unnest() %>% ungroup()
}else if(Big == 'B'){
# =---------------------------------------------------------
data_split <- data_split %>% filter(time == 'future') %>%
group_by(id, ISO3, county, Country, x, y, time, Big, gSeason) %>%
summarise_all(mean) %>%
mutate_at(.vars = vars(SLGP, LGP), .funs = function(x){round(x, 0)}) %>%
ungroup() %>%
bind_rows(filter(data_split, time == 'past'), .)
# gSeason...
gSeason_i <- data_split %>% dplyr::select( time, x, y, gSeason) %>%
group_by( time, x, y) %>% summarise(gSeason = max(gSeason, na.rm = TRUE)) %>%
ungroup() %>% group_by(time, x, y) %>%
summarise(gSeason = round(mean(gSeason, na.rm = TRUE), 0)) %>% ungroup()
# =----------------------------------------------------------------------
# gSeason...
two_index <- data_split %>% dplyr::select( time, x, y, gSeason, SLGP, LGP) %>%
group_by(time, x, y, gSeason) %>%
summarise_all(.f = function(x){round(mean(x, na.rm = TRUE), 0)}) %>%
# ungroup(year) %>%
arrange(gSeason) %>%
tidyr::nest(data = c('x', 'y','SLGP', 'LGP')) %>%
drop_na() %>% filter(gSeason < 3) %>% tidyr::unnest() %>% ungroup()
} else{print('Change big argument... >.<')}
# =---------------------------------------------------------------------
# gSeason
limits_gs <- dplyr::select(gSeason_i, gSeason) %>% summarise_all(.funs = c('min', 'max'))
# # Primero dejar? hechos los de presente... luego repito los de futuro...
gs <- ggplot() +
geom_tile(data = filter(gSeason_i, time == 'past'), aes(x = x, y = y, fill = gSeason)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('gSeason\n(day) '), title = 'Historic', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = as.integer(limits_gs),
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/gSeason_past.png') , width = 10, height = 10, dpi = 300)
# # =- Futuro.
gs_f <- ggplot() +
geom_tile(data = filter(gSeason_i, time == 'future'), aes(x = x, y = y, fill = gSeason)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('gSeason\n(day) '), title = 'Future', x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = as.integer(limits_gs),
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/gSeason_future.png') , width = 10, height = 10, dpi = 300)
# =--------
gSeason_dif <- gSeason_i %>% pivot_wider(names_from = time, values_from = gSeason) %>%
mutate(gSeason = future - past)
c <- ggplot() +
geom_tile(data = gSeason_dif, aes(x = x, y = y, fill = gSeason)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('gSeason\n(days) '), title = 'Change', x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(low = '#A50026', mid = 'white', high = '#000099',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_gSeason.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/all_gSeason.png'), width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(gs, gs_f, c, ncol=3))
dev.off()
# =---------------------------------------------------------------------
# =---------------------------------------------------------------------
# Aqui hay que hacer los filtros de acuerdo al grupo ...
# rows_n <- data_split %>% filter(time == 'past', year == '1985') %>% dplyr::select(x, y, id) %>% unique() %>% nrow()
SLGP_dif <- two_index %>% dplyr::select(-LGP) %>%
pivot_wider(names_from = time, values_from = SLGP) %>%
mutate(SLGP = future - past, gSeason = glue::glue('gSeason = {gSeason}'))
limits_two <- two_index %>% group_by(gSeason) %>%
dplyr::select(SLGP, LGP) %>% summarise_all(.funs = c('min', 'max')) %>%
ungroup()
gS <- unique(two_index$gSeason)
for(i in gS){
my_limits <- c(limits_two$SLGP_min[i], limits_two$SLGP_max[i])
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
SLGP_p_1 <- ggplot(filter(two_index %>% mutate(gSeason = glue::glue('gSeason = {gSeason}')),
time == 'past', gSeason == glue::glue('gSeason = {i}'))) +
geom_tile(aes(x = x, y = y, fill = SLGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('SLGP\n(Day of\nthe year) '),
title = glue::glue('gSeason = {i}; Historic'),
x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/SLGP_past_{i}.png') , width = 10, height = 10, dpi = 300)
SLGP_f_1 <- ggplot(filter(two_index %>% mutate(gSeason = glue::glue('gSeason = {gSeason}')),
time == 'future', gSeason == glue::glue('gSeason = {i}'))) +
geom_tile(aes(x = x, y = y, fill = SLGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('SLGP\n(Day of\nthe year) '),
title = glue::glue('gSeason = {i}; Future'), x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/SLGP_future_{i}.png') , width = 10, height = 10, dpi = 300)
# =--------
dt_MM <- filter(SLGP_dif, gSeason == glue::glue('gSeason = {i}'))$SLGP
my_limits_MM <- c(min(dt_MM), max(dt_MM))
my_breaks_MM <- round(seq(my_limits_MM[1], my_limits_MM[2], length.out= 3), 0)
my_limits_MM <- c(ifelse(my_limits_MM[1] > my_breaks_MM[1], my_breaks_MM[1], my_limits_MM[1]) ,ifelse(my_limits_MM[2] < my_breaks_MM[3], my_breaks_MM[3], my_limits_MM[2]))
d <- ggplot() +
geom_tile(data = filter(SLGP_dif, gSeason == glue::glue('gSeason = {i}')),
aes(x = x, y = y, fill = SLGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('SLGP\n(days) '),
title = glue::glue('gSeason = {i}; Change'),x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2(limits = my_limits_MM,
breaks = my_breaks_MM,
low = '#A50026', mid = 'white', high = '#000099',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_SLGP_{i}.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/all_SLGP_{i}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(SLGP_p_1, SLGP_f_1, d, ncol=3))
dev.off()
}
#### =---------------------------------------------------------------------------
# =--------------------------------------------------------------------
# =--------------------------------------------------------------------
LGP_dif <- two_index %>% dplyr::select(-SLGP) %>%
pivot_wider(names_from = time, values_from = LGP) %>%
mutate(LGP = future - past, gSeason = glue::glue('gSeason = {gSeason}'))
for(i in gS){
my_limits <- c(limits_two$LGP_min[i], limits_two$LGP_max[i])
my_breaks <- round(seq(my_limits[1], my_limits[2], length.out= 3), 0)
my_limits <- c(ifelse(my_limits[1] > my_breaks[1], my_breaks[1], my_limits[1]) ,ifelse(my_limits[2] < my_breaks[3], my_breaks[3], my_limits[2]))
LGP_p <- ggplot(filter(two_index %>% mutate(gSeason = glue::glue('gSeason = {gSeason}')),
time == 'past', gSeason == glue::glue('gSeason = {i}'))) +
geom_tile(aes(x = x, y = y, fill = LGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('LGP\n(days) '),
title = glue::glue('gSeason = {i}; Historic'),x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/LGP_past_{i}.png') , width = 10, height = 10, dpi = 300)
LGP_f <- ggplot(filter(two_index %>% mutate(gSeason = glue::glue('gSeason = {gSeason}')),
time == 'future', gSeason == glue::glue('gSeason = {i}'))) +
geom_tile(aes(x = x, y = y, fill = LGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('LGP\n(days) '),
title = glue::glue('gSeason = {i}; Future'), x = 'Longitude', y = 'Latitude') +
scale_fill_viridis_c(limits = my_limits,
breaks = my_breaks,
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/LGP_future_{i}.png') , width = 10, height = 10, dpi = 300)
# =--------
dt_MM <- filter(LGP_dif, gSeason == glue::glue('gSeason = {i}'))$LGP
my_limits_MM <- c(min(dt_MM), max(dt_MM))
my_breaks_MM <- round(seq(my_limits_MM[1], my_limits_MM[2], length.out= 3), 0)
my_limits_MM <- c(ifelse(my_limits_MM[1] > my_breaks_MM[1], my_breaks_MM[1], my_limits_MM[1]) ,ifelse(my_limits_MM[2] < my_breaks_MM[3], my_breaks_MM[3], my_limits_MM[2]))
e <- ggplot(filter(LGP_dif, gSeason == glue::glue('gSeason = {i}')) ) +
geom_tile(aes(x = x, y = y, fill = LGP)) +
geom_sf(data = glwd1, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = glwd2, fill = 'lightblue', color = 'lightblue') +
geom_sf(data = pais, fill = NA, color = gray(.8)) +
geom_sf(data = shp_sf, fill = NA, color = gray(.1)) +
coord_sf(xlim = xlims, ylim = ylims) +
labs(fill = glue::glue('LGP\n(days) '),
title = glue::glue('gSeason = {i}; Change'),x = 'Longitude', y = 'Latitude') +
scale_fill_gradient2( limits = my_limits_MM,
breaks = my_breaks_MM,
low = '#A50026', mid = 'white', high = '#000099',
guide = guide_colourbar(barwidth = 20, label.theme = element_text(angle = 25, size = 35))) +
scale_y_continuous(breaks = round(ylims, 2), n.breaks = 3) +
scale_x_continuous(breaks = round(xlims, 2), n.breaks = 3) +
theme_bw() + theme(legend.position = 'bottom', text = element_text(size=35),
legend.title=element_text(size=35),
legend.spacing = unit(5, units = 'cm'),
legend.spacing.x = unit(1.0, 'cm'), plot.title = element_text(hjust = 0.5))
ggsave(glue::glue('{path}{Country}/graphs/{county}/maps/Dif_LGP_{i}.png') , width = 10, height = 10, dpi = 300)
png(filename = glue::glue('{path}{Country}/graphs/{county}/maps/all_LGP_{i}.png') , width=28,height=10,units="in", res = 300)
print(gridExtra::grid.arrange(LGP_p, LGP_f, e, ncol=3))
dev.off()
}
}
## =-----------------------------------------
Clim_graph <- function(historic){
Country <- unique(historic$Country)
shp_sf <- shp %>% sf::st_as_sf()
country <- country1 %>% sf::st_as_sf()
xlims <- sf::st_bbox(shp_sf)[c(1, 3)]
ylims <- sf::st_bbox(shp_sf)[c(2, 4)]
# # =------------------------------------------
# glwd1 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_1.shp' )
# crs(glwd1) <- crs(shp)
# glwd1 <- rgeos::gSimplify(glwd1, tol = 0.05, topologyPreserve = TRUE) %>%
# sf::st_as_sf()
#
# glwd2 <- raster::shapefile('//dapadfs/workspace_cluster_8/climateriskprofiles/data/shps/GLWD/glwd_2.shp' )
# crs(glwd2) <- crs(shp)
# glwd2 <- rgeos::gSimplify(glwd2, tol = 0.05, topologyPreserve = TRUE) %>%
# sf::st_as_sf()