rm(list=ls())
library(Hmisc)
library(data.table)
library(geosphere)

setwd('~/research/gTiger/Article_citsci/open_code')
source('settings/parameters.r')
source('settings/paths.r')
source('settings/functions.r')

setwd(path_to_data)

D = readRDS('mosquito_alert_model_dataset.Rds')

D = D[SE>0]

class(D$biweek_start_date)

D$masked_lon = as.numeric(sapply(D$TigacellID, function(x) strsplit(x, "_")[[1]][1]))
D$masked_lat = as.numeric(sapply(D$TigacellID, function(x) strsplit(x, "_")[[1]][2]))

D$diag_dist_m = distVincentyEllipsoid(D[, list(masked_lon, masked_lat)], D[, list(masked_lon, masked_lat)] + .05)
D$area_sq_km = ((D$diag_dist_m^2)/2)/(1000*1000)

MA_biweek_sq_kms = sum(D$area_sq_km)

MA_cost_per_biweek_sq_km = 300000/MA_biweek_sq_kms


MA_cost_per_month_sq_km = 2*MA_cost_per_biweek_sq_km

D = readRDS('ovitrap_model_dataset.Rds')

D = D[NAME_2=="Murcia"]

D$biweek = as.integer(difftime(D$F_date, season_day_0_2014, units='weeks')/2)

D$Tigacell_biweek = paste0(D$TigacellID, "_", D$biweek)

D$masked_lon = as.numeric(sapply(D$TigacellID, function(x) strsplit(x, "_")[[1]][1]))
D$masked_lat = as.numeric(sapply(D$TigacellID, function(x) strsplit(x, "_")[[1]][2]))

D$diag_dist_m = distVincentyEllipsoid(D[, list(masked_lon, masked_lat)], D[, list(masked_lon, masked_lat)] + .05)
D$area_sq_km = ((D$diag_dist_m^2)/2)/(1000*1000)

D_agg = D[ , list(area_sq_km=last(area_sq_km)), by=Tigacell_biweek]

ovi_biweek_sq_km = sum(D_agg$area_sq_km)

ovi_cost_per_biweek_sq_km = 2*90537/ovi_biweek_sq_km

ovi_cost_per_month_sq_km = 2*ovi_cost_per_biweek_sq_km



ovi_cost_per_month_sq_km/ MA_cost_per_month_sq_km 


MA_cost_per_month_sq_km /ovi_cost_per_month_sq_km
######


make_trap_map = function(){
par(mar=c(0.1,0.1,0.1,0.1))
plot(gadm_prov_utm, bg="#000066",col="#000033")
plot(D_utm, pch=20, col="#ffff9977", cex=.1, add=T)
}

png('Article_citsci/article/figures/ovitrap_locations_2014-15.png')
make_trap_map()
dev.off()

pdf('Article_citsci/article/figures/ovitrap_locations_2014-15.pdf')
make_trap_map()
dev.off()


 (64368+5040+2850)/74824 

