library(tidyverse)
library(scales)
library(janitor)
library(paletteer)
source("Scripts/Return_Data_Center_Model_Inputs.r")
DATA_CENTER <- DATA_CENTER_INPUTS()
DATA_CENTER <- rbind(DATA_CENTER,do.call(rbind, lapply(2031:2050,function(x){DATA_CENTER  %>% filter(Year==2030) %>% mutate(Year=x)})))
GET_DATA_CENTER_RES <- function(DIR,EVENT,REGION='Wyoming'){
#	DIR<-"Data_Center_Employment_100_People"
#	EVENT<-'DATA_CENTER_EMP'
	PATH <- paste0("Model_Outputs/IMPLAN/Data_Center/",REGION,"/",DIR,"/economic_indicators_by_impact.csv")
	INPUT_DATA <- read_csv(PATH) %>% rename(Income=`Labor Income`)  %>% mutate(Income=parse_number(Income),`Value Added`=parse_number(`Value Added`),Output=parse_number(Output))  %>% mutate(Region=REGION,Event=EVENT)
	INPUT_DATA$Impact <- gsub("3 - ","",gsub("2 - ","",gsub("1 - ","",INPUT_DATA$Impact)))
	INPUT_DATA <-  INPUT_DATA %>% pivot_longer(-c(Impact,Region,Event)) %>% rename(Type=Impact,Impact=name)
	INPUT_DATA$Impact <- gsub("3 - ","",gsub("2 - ","",gsub("1 - ","",INPUT_DATA$Impact)))
	INPUT_DATA <- INPUT_DATA %>% mutate(Source='IMPLAN') %>% select(Type,Region,Impact,Event,Source,value)
	INPUT_DATA <- INPUT_DATA %>% filter(!is.na(Type))
	INPUT_DATA$value <- if(EVENT=='DATA_CENTER_EMP'){INPUT_DATA$value <- INPUT_DATA$value/100}else{INPUT_DATA$value/10^6} #shift to impact per dollar from impact per million dollar
INPUT_DATA <-	INPUT_DATA %>% rename('marginal'=value)
INPUT_DATA <- INPUT_DATA   %>% left_join(DATA_CENTER) %>% mutate(value=marginal*total)
INPUT_DATA <- INPUT_DATA %>% group_by(Type,Year,Region,Impact,Source,Event) %>% summarize(value=sum(value)) %>% ungroup
return(INPUT_DATA )
}

WY_DATA <- rbind(GET_DATA_CENTER_RES("Computer_Purchase_Mill_Dollars",'COMPUTERS'),GET_DATA_CENTER_RES("Construction_Mill_Dollars",'BUILDING'),GET_DATA_CENTER_RES("Data_Center_Employment_100_People",'DATA_CENTER_EMP'),GET_DATA_CENTER_RES("Other_Purchases_Mill_Dollars","OTHER"),GET_DATA_CENTER_RES("Pre-Development_Million_Dollars","PRE_DEV"),GET_DATA_CENTER_RES("Utlities_Electricity_Million_USD","UTILITIES"))
US_DATA <- rbind(GET_DATA_CENTER_RES("Computer_Purchase_Mill_Dollars",'COMPUTERS','US'),GET_DATA_CENTER_RES("Construction_Mill_Dollars",'BUILDING','US'),GET_DATA_CENTER_RES("Data_Center_Employment_100_People",'DATA_CENTER_EMP','US'),GET_DATA_CENTER_RES("Other_Purchases_Mill_Dollars","OTHER",'US'),GET_DATA_CENTER_RES("Pre-Development_Million_Dollars","PRE_DEV",'US'),GET_DATA_CENTER_RES("Utlities_Electricity_Million_USD","UTILITIES",'US'))
DATA <-rbind(US_DATA,WY_DATA) %>% group_by(Year,Type,Region,Impact,Source) %>% summarize(value=sum(value)) %>% ungroup

######################
GET_DATA_CENTER_TAX <- function(DIR,EVENT,REGION='Wyoming'){
#	DIR<-"Data_Center_Employment_100_People"
#	EVENT<-'DATA_CENTER_EMP'
	PATH <- paste0("Model_Outputs/IMPLAN/Data_Center/",REGION,"/",DIR,"/tax_results.csv")
	INPUT_DATA <- read_csv(PATH)  %>% mutate(`Sub County General`=parse_number(`Sub County General`),`Sub County Special Districts`=parse_number(`Sub County Special Districts`),County=parse_number(County),County=County+`Sub County General`+`Sub County Special Districts`,State=parse_number(State),Federal=parse_number(Federal)) %>% select(-Total,`Sub County General`,`Sub County Special Districts`)  %>% mutate(Region=REGION,Event=EVENT)
	INPUT_DATA$Impact <- gsub("3 - ","",gsub("2 - ","",gsub("1 - ","",INPUT_DATA$Impact)))
	INPUT_DATA <-  INPUT_DATA %>% pivot_longer(-c(Impact,Region,Event)) %>% rename(Type=Impact,Impact=name)
	INPUT_DATA$Impact <- gsub("3 - ","",gsub("2 - ","",gsub("1 - ","",INPUT_DATA$Impact)))
	INPUT_DATA <- INPUT_DATA %>% mutate(Source='IMPLAN') %>% select(Type,Region,Impact,Event,Source,value)
	INPUT_DATA <- INPUT_DATA %>% filter(!is.na(Type))
	INPUT_DATA$value <- if(EVENT=='DATA_CENTER_EMP'){INPUT_DATA$value <- INPUT_DATA$value/100}else{INPUT_DATA$value/10^6} #shift to impact per dollar from impact per million dollar
INPUT_DATA <-	INPUT_DATA %>% rename('marginal'=value)
INPUT_DATA <- INPUT_DATA   %>% left_join(DATA_CENTER) %>% mutate(value=marginal*total)
INPUT_DATA <- INPUT_DATA %>% group_by(Type,Year,Region,Impact,Source,Event) %>% summarize(value=sum(value)) %>% ungroup
return(INPUT_DATA )
}
######################
#REGION='Wyoming'
#DIR <- "Computer_Purchase_Mill_Dollars"
#EVENT <- 'COMPUTERS'
GET_DATA_CENTER_JOBS <- function(DIR,EVENT,REGION='Wyoming'){
#	DIR<-"Data_Center_Employment_100_People"
#	EVENT<-'DATA_CENTER_EMP'
	PATH <- paste0("Model_Outputs/IMPLAN/Data_Center/",REGION,"/",DIR,"/occupation_impacts_table.csv")
	INPUT_DATA <- read_csv(PATH)[,c(2:3,5)]
	colnames(INPUT_DATA) <- c("Occupation","Employment","Wages")
	INPUT_DATA$Wages <- parse_number(INPUT_DATA$Wages)

	INPUT_DATA$Employment<- if(EVENT=='DATA_CENTER_EMP'){INPUT_DATA$Employment<- INPUT_DATA$Employment/100}else{INPUT_DATA$Employment/10^6} #shift to impact per dollar from impact per million dollar

	INPUT_DATA$Wages<- if(EVENT=='DATA_CENTER_EMP'){INPUT_DATA$Wages<- INPUT_DATA$Wages/100}else{INPUT_DATA$Wages/10^6} #shift to impact per dollar from impact per million dollar
INPUT_DATA <-	INPUT_DATA %>% rename('marginal_Employment'=Employment,'marginal_Wages'=Wages)
INPUT_DATA$Event <-  EVENT
INPUT_DATA$Region <- REGION 
INPUT_DATA <- 	INPUT_DATA   %>% left_join(DATA_CENTER) %>% mutate(Employment=marginal_Employment*total,Wages=marginal_Wages*total)
INPUT_DATA <- INPUT_DATA %>% group_by(Occupation,Year,Region,Event) %>% summarize(Employment=sum(Employment),Wages=sum(Wages)) %>% ungroup
return(INPUT_DATA )
}

WY_DATA_TAX <- rbind(GET_DATA_CENTER_TAX("Computer_Purchase_Mill_Dollars",'COMPUTERS'),GET_DATA_CENTER_TAX("Construction_Mill_Dollars",'BUILDING'),GET_DATA_CENTER_TAX("Data_Center_Employment_100_People",'DATA_CENTER_EMP'),GET_DATA_CENTER_TAX("Other_Purchases_Mill_Dollars","OTHER"),GET_DATA_CENTER_TAX("Pre-Development_Million_Dollars","PRE_DEV"),GET_DATA_CENTER_TAX("Utlities_Electricity_Million_USD","UTILITIES"))
US_DATA_TAX <- rbind(GET_DATA_CENTER_TAX("Computer_Purchase_Mill_Dollars",'COMPUTERS','US'),GET_DATA_CENTER_TAX("Construction_Mill_Dollars",'BUILDING','US'),GET_DATA_CENTER_TAX("Data_Center_Employment_100_People",'DATA_CENTER_EMP','US'),GET_DATA_CENTER_TAX("Other_Purchases_Mill_Dollars","OTHER",'US'),GET_DATA_CENTER_TAX("Pre-Development_Million_Dollars","PRE_DEV",'US'),GET_DATA_CENTER_TAX("Utlities_Electricity_Million_USD","UTILITIES",'US'))
TAX <-rbind(US_DATA_TAX,WY_DATA_TAX) %>% group_by(Year,Type,Region,Impact,Source) %>% summarize(value=sum(value)) %>% ungroup
IMPLAN <- rbind(TAX,DATA ) %>% select(Type,Year,Region,Impact,Source,value)
IMPLAN$Region <- ifelse(IMPLAN$Region=='Wyoming','WY',IMPLAN$Region)
########################
IMPLAN_DOLLAR <- IMPLAN %>% filter(Impact!='Employment')
IMPLAN_DOLLAR <-  IMPLAN_DOLLAR  %>% group_by(Type,Year,Impact,Source) %>% mutate(value=ifelse(Region=='US',value-min(value),value) ) %>% ungroup
#State tax write off
IMPLAN_DOLLAR[IMPLAN_DOLLAR$Impact %in% c("County","State") & IMPLAN_DOLLAR$Type=='Direct'& IMPLAN_DOLLAR$Region=='WY',"value"] <- 0


#TEMP <- IMPLAN_DOLLAR %>% filter(Year==2028)  %>% mutate(value=0)
#IMPLAN_DOLLAR <- rbind(IMPLAN_DOLLAR,rbind(TEMP %>% mutate(Year=2025),TEMP %>% mutate(Year=2026),TEMP %>% mutate(Year=2027)))

	TAX_STATE_BASE <- 500*(10^6*0.115*((62.835-12)/1000)) #State
	TAX_COUNTY_BASE <-500*(10^6*0.115*((12)/1000)) #County
	IMPLAN_DOLLAR[IMPLAN_DOLLAR$Impact %in% c("County") & IMPLAN_DOLLAR$Type=='Direct'& IMPLAN_DOLLAR$Region=='WY',]
	DIRECT_TAX_FLOW <- (c(seq(0,1,by=1/3),rev(seq(0,1,by=1/40))))[1:25]
	IMPLAN_DOLLAR[IMPLAN_DOLLAR$Impact %in% c("County") & IMPLAN_DOLLAR$Type=='Direct'& IMPLAN_DOLLAR$Region=='WY','value']	<- DIRECT_TAX_FLOW *TAX_COUNTY_BASE 
	IMPLAN_DOLLAR[IMPLAN_DOLLAR$Impact %in% c("State") & IMPLAN_DOLLAR$Type=='Direct'& IMPLAN_DOLLAR$Region=='WY','value'] <- DIRECT_TAX_FLOW *TAX_STATE_BASE 


IMPLAN_DOLLAR <- 	IMPLAN_DOLLAR %>% filter(!(Impact %in% c("Sub County Special Districts","Sub County County General")))


IMPLAN_DOLLAR$Region <- factor(ifelse(IMPLAN_DOLLAR$Region=="WY",'Wyoming',"Rest of United States"),levels=rev(c("Wyoming","Rest of United States")))
IMPLAN_DOLLAR_REGION <- IMPLAN_DOLLAR %>% filter(Impact %in% c('Income','Output','Value Added')) %>% group_by(Year,Region,Impact) %>% summarize(value=sum(value)) %>% ungroup

########################Plot of total impact in US and in Wyoming
IMPLAN_REGION_DOLLAR_PLOT <- ggplot(IMPLAN_DOLLAR_REGION ,aes(x=Year,y=value/10^9,fill=Region,group=Region))+geom_area()+facet_wrap(~Impact,ncol=1)+ylab("Billion Dollars")+theme_bw()+ theme(plot.title = element_text(size = 22, face = "bold"), axis.title = element_text(size = 18),axis.text  = element_text(size = 16),legend.text = element_text(size = 14),strip.text = element_text(size = 14),legend.position="top")+  theme(legend.title=element_blank())+scale_fill_manual(values=c("#3C3B6E","#FFC425"))+scale_x_continuous(breaks=seq(2025,2060,by=5))+scale_y_continuous(breaks=seq(0,5000,by=1))
ggsave( filename = "./Results/Data Center: Dollar Values of Model.png",  plot = IMPLAN_REGION_DOLLAR_PLOT,  width = 8.5,  height = 11*(3./4),  units = "in", dpi = 300)
########################Taxes
IMPLAN_TAX <- IMPLAN_DOLLAR %>% filter(Impact %in% c("County","State","Federal"))
OLD <- IMPLAN_TAX 
IMPLAN_TAX <- IMPLAN_TAX %>% group_by(Type,Year,Impact,Source) %>% mutate(value=ifelse(Region=='Rest of United States',value-min(value),value) ) %>% ungroup

IMPLAN_TAX <- IMPLAN_TAX %>% mutate(Impact=ifelse(Region=='Rest of United States' & Impact!='Federal',"Other States",as.character(Impact)))
IMPLAN_TAX$Impact <- ifelse(IMPLAN_TAX$Impact =='State',"Wyoming State Taxes",IMPLAN_TAX$Impact)
IMPLAN_TAX$Impact <- ifelse(IMPLAN_TAX$Impact =='County',"Wyoming County Taxes",IMPLAN_TAX$Impact)
IMPLAN_TAX$Impact <- factor(IMPLAN_TAX$Impact,levels=c('Output','Value Added','Income','Federal','Other States','Wyoming State Taxes','Wyoming County Taxes'))

#IMPLAN_TAX <- IMPLAN_TAX[!(IMPLAN_TAX$Type=='Direct' & IMPLAN_TAX$Region=='Rest of United States'& IMPLAN_TAX$Impact!='Federal'),] 
#IMPLAN_TAX <-IMPLAN_TAX[!(IMPLAN_TAX$Type=='Direct' & IMPLAN_TAX$Region=='Wyoming'& IMPLAN_TAX$Impact=='Federal'),] 
IMPLAN_DOLLAR <- IMPLAN_DOLLAR %>% filter(!(Impact %in% c("County","State","Federal","Other States"))) %>% rbind(IMPLAN_TAX)
IMPLAN_TAX_FIG_DATA <- IMPLAN_TAX  %>% group_by(Year,Impact) %>% summarize(value=sum(value)) %>% ungroup
IMPLAN_TAX %>% filter(Impact=='Other States',Type=='Direct',value>0)
TAX_PLOT <- ggplot(IMPLAN_TAX_FIG_DATA, aes(x=Year,y=value/10^6,fill=Impact,group=Impact))+geom_area()+ylab("Million Dollars")+theme_bw()+ theme(plot.title = element_text(size = 22, face = "bold"), axis.title = element_text(size = 18),axis.text  = element_text(size = 16),legend.text = element_text(size = 14),strip.text = element_text(size = 14),legend.position="top")+  theme(legend.title=element_blank())+scale_fill_manual(values=c("#3C3B6E","firebrick","#FFC425","#492F24"))+scale_x_continuous(breaks=seq(2025,2060,by=5))+scale_y_continuous(breaks=seq(0,5000,by=50))


#############Dollar values table
NPV <- IMPLAN_DOLLAR %>% mutate(value=value/10^9,Discount_0=value,Discount_2=value/(1+0.02)^(Year-2026),Discount_5=value/(1+0.05)^(Year-2026),Discount_10=value/(1+0.1)^(Year-2026)) %>% group_by(Impact,Region) %>% summarize(Discount_0=sum(Discount_0),Discount_2=sum(Discount_2),Discount_5=sum(Discount_5),Discount_10=sum(Discount_10)) %>% ungroup
NPV[,3:6] <-round(NPV[,3:6],2)
NPV
WY_NPV <- cbind(c("Discount","0%","2%","5%","10%"),NPV %>% filter(Region=='Wyoming') %>% select(-Region) %>% t) %>% as_tibble
colnames(WY_NPV ) <- WY_NPV[1,]
WY_NPV <- WY_NPV[-1,] 
WY_NPV <- WY_NPV %>%  mutate_if(is.character, str_trim)
WY_NPV[,-1] <- lapply(	lapply(WY_NPV[,-1] ,as.numeric) %>% as_tibble,dollar) %>% as_tibble
WY_TAX <- WY_NPV[,c(1,5:7)]
WY_NPV <- WY_NPV[,c(1:4)]


#US
US_NPV <- cbind(c("Discount","0%","2%","5%","10%"),NPV %>% filter(Region!='Wyoming') %>% select(-Region) %>% t) %>% as_tibble
colnames(US_NPV ) <- US_NPV[1,]
US_NPV <- US_NPV[-1,] 
US_NPV <- US_NPV %>%  mutate_if(is.character, str_trim)
US_NPV[,-1] <- lapply(	lapply(US_NPV[,-1] ,as.numeric) %>% as_tibble,dollar) %>% as_tibble
US_TAX <- US_NPV[,c(1,5:6)]
US_NPV <- US_NPV[,c(1:4)]
ALL_TAX <- US_TAX
ALL_TAX[,2] <- dollar(parse_number(as.character(US_TAX[,2] %>% t))+ parse_number(as.character(WY_TAX[,2] %>% t)))
ALL_TAX <- cbind(ALL_TAX,WY_TAX[,-1:-2])

#write.csv(US_TAX,"./Results/US_Tax_Values.csv",row.names=FALSE)



############Employment
IMPLAN_EMPLOY <- IMPLAN %>% filter(Impact=='Employment',Source=='IMPLAN') %>% group_by(Year,Region,Impact) %>% summarize(value=sum(value)) %>% ungroup
IMPLAN_EMPLOY <-  IMPLAN_EMPLOY  %>% group_by(Year,Impact) %>% mutate(value=ifelse(Region=='US',value-min(value),value) ) %>% ungroup
##########
	EMP_MAIN_SUMMARY <- IMPLAN %>% filter(Impact=='Employment',Source=='IMPLAN')  %>% group_by(Year,Region,Type)%>% summarize(value=mean(value)) %>% group_by(Year,Type)
US <- EMP_MAIN_SUMMARY  %>% filter(Region=='US') %>% rename(US_value=value) %>% ungroup %>% left_join(EMP_MAIN_SUMMARY  %>% filter(Region=='WY') %>% select(-Region)) %>% mutate(US_value=US_value-value,US_value=ifelse(US_value<0,0,US_value)) %>% select(-value) %>% rename(value=US_value)
EMP_MAIN_SUMMARY <-  rbind(US,EMP_MAIN_SUMMARY  %>% filter(Region=='WY'))%>% mutate(Period=ifelse(Year>CONST_END_YEAR,'Operations',"Construction")) %>% group_by(Period,Region,Type) %>% summarize(value=mean(value)) %>% ungroup %>% pivot_wider(values_from=value,names_from=c(Period)) %>% arrange(desc(Region))

EMP_MAIN_SUMMARY$Region <- ifelse(EMP_MAIN_SUMMARY$Region=='WY','Wyoming','Other States')
EMP_MAIN_SUMMARY <- rbind(EMP_MAIN_SUMMARY,c("United States","Total",as.numeric(colSums(EMP_MAIN_SUMMARY[-1:-2]))))
EMP_MAIN_SUMMARY[,3] <- round(as.numeric(t(EMP_MAIN_SUMMARY[,3] )),0)
EMP_MAIN_SUMMARY[,4] <- round(as.numeric(t(EMP_MAIN_SUMMARY[,4] )),0)


##################
IMPLAN_EMPLOY$Region <- factor(ifelse(IMPLAN_EMPLOY$Region=="WY",'Wyoming',"Rest of United States"),levels=rev(c("Wyoming","Rest of United States")))
 IMPLAN_EMPLOYMENT_PLOT <- ggplot(IMPLAN_EMPLOY ,aes(x=Year,y=value,fill=Region,group=Region))+geom_area()+ylab("Employment (Full-Time Equivalent)")+theme_bw()+ theme(plot.title = element_text(size = 22, face = "bold"), axis.title = element_text(size = 18),axis.text  = element_text(size = 16),legend.text = element_text(size = 14),strip.text = element_text(size = 14),legend.position="top")+  theme(legend.title=element_blank())+scale_fill_manual(values=c("#3C3B6E","#FFC425"))+scale_x_continuous(breaks=seq(2025,2060,by=5))+scale_y_continuous(breaks=c(seq(0,50000,by=1000)),labels = label_comma())
 
####Occupation
CLEAN_NAMES <- function(DATASET){
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Doctors','1Doctors',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Maintenance workers','1Maintenance workers',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Restaurants','1Restaurants',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Business Operations','1Business Operations',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Metal Workers','1Metal Workers',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Material Moving','1Material Moving',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='General Production','1General Production',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Management','1Management',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='Clerks','1Clerks',DATASET$Occupation)
#DATASET$Occupation <- ifelse(DATASET$Occupation=='System Operators','1System Operators',DATASET$Occupation)

DATASET$Occupation <- ifelse(DATASET$Occupation=='Healthcare Diagnosing or Treating Practitioners','Doctors',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Other Installation, Maintenance, and Repair Occupations','Maintenance workers',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Food and Beverage Serving Workers','Restaurants',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Business Operations Specialists','Business Operations',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Metal Workers and Plastic Workers','Metal Workers',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Material Moving Workers','Material Moving',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Other Production Occupations','General Production',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Other Management Occupations','Management',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Information and Record Clerks','Clerks',DATASET$Occupation)
DATASET$Occupation <- ifelse(DATASET$Occupation=='Plant and System Operators','System Operators',DATASET$Occupation)
return(DATASET)
}

WY_DATA_OCCUPATION <- rbind(GET_DATA_CENTER_JOBS("Computer_Purchase_Mill_Dollars",'COMPUTERS'),GET_DATA_CENTER_JOBS("Construction_Mill_Dollars",'BUILDING'),GET_DATA_CENTER_JOBS("Data_Center_Employment_100_People",'DATA_CENTER_EMP'),GET_DATA_CENTER_JOBS("Other_Purchases_Mill_Dollars","OTHER"),GET_DATA_CENTER_JOBS("Pre-Development_Million_Dollars","PRE_DEV"),GET_DATA_CENTER_JOBS("Utlities_Electricity_Million_USD","UTILITIES"))
WY_DATA_OCCUPATION <- CLEAN_NAMES(WY_DATA_OCCUPATION )

US_DATA_OCCUPATION <- rbind(GET_DATA_CENTER_JOBS("Computer_Purchase_Mill_Dollars",'COMPUTERS','US'),GET_DATA_CENTER_JOBS("Construction_Mill_Dollars",'BUILDING','US'),GET_DATA_CENTER_JOBS("Data_Center_Employment_100_People",'DATA_CENTER_EMP','US'),GET_DATA_CENTER_JOBS("Other_Purchases_Mill_Dollars","OTHER",'US'),GET_DATA_CENTER_JOBS("Pre-Development_Million_Dollars","PRE_DEV",'US'),GET_DATA_CENTER_JOBS("Utlities_Electricity_Million_USD","UTILITIES",'US'))
US_DATA_OCCUPATION <- CLEAN_NAMES(US_DATA_OCCUPATION )




RES <- rbind(US_DATA_OCCUPATION ,WY_DATA_OCCUPATION) %>% group_by(Occupation,Year,Region) %>% summarize(Employment=sum(Employment),Wages=sum(Wages)) %>% ungroup

JOBS <- rbind(RES %>% filter(Year<2030,Region=='Wyoming') %>% group_by(Occupation) %>% summarize(Peak_Employment=max(Employment,na.rm=TRUE),Employment=mean(Employment,na.rm=TRUE),Wages=mean(Wages,na.rm=TRUE)) %>% ungroup %>% arrange(desc(Employment)) %>% mutate(Region='Wyoming',Period='Construction'),RES %>% filter(Year>=2030,Region=='Wyoming') %>% arrange(desc(Employment)) %>% mutate(Peak_Employment=Employment,Region='Wyoming',Period='Operations')%>% select(-Year),RES %>% filter(Year<2030,Region!='Wyoming') %>% group_by(Occupation) %>% summarize(Peak_Employment=max(Employment,na.rm=TRUE),Employment=mean(Employment,na.rm=TRUE),Wages=mean(Wages,na.rm=TRUE)) %>% ungroup %>% arrange(desc(Employment)) %>% mutate(Region='United States',Period='Construction'),RES %>% filter(Year>=2030,Region!='Wyoming') %>% arrange(desc(Employment)) %>% mutate(Peak_Employment=Employment,Region='United States',Period='Operations')%>% select(-Year))

US_DATA_OCCUPATION <- CLEAN_NAMES(US_DATA_OCCUPATION )

RES <- WY_DATA_OCCUPATION %>% select(Occupation,Employment,Wages,Year)
RES_US <- US_DATA_OCCUPATION  %>% select(Occupation,Employment,Wages,Year)




###
WY_JOBS <- JOBS %>% filter(Region=='Wyoming') %>%  rename(Peak_Wy_Employment=Peak_Employment,Wy_Employment=Employment,Wy_Wages=Wages)  %>% select(-Region)
US_JOBS <- JOBS %>% filter(Region!='Wyoming') %>%  rename(Peak_US_Employment=Peak_Employment,US_Employment=Employment,US_Wages=Wages)  %>% select(-Region)

JOINED <- WY_JOBS %>% left_join(US_JOBS)  %>% mutate(US_Employment=US_Employment-Wy_Employment,Peak_US_Employment=Peak_Wy_Employment,US_Wages=US_Wages-Wy_Wages) 
JOINED$US_Employment  <-  ifelse(JOINED$US_Employment<0,0,JOINED$US_Employment)
JOINED$US_Wages <-  ifelse(JOINED$US_Employment==0,0,JOINED$US_Wages)
JOINED$Peak_US_Employment <-  ifelse(JOINED$US_Employment==0,0,JOINED$Peak_US_Employment)

EMP_SUMMARY <- JOINED %>% mutate(Total_Employment=US_Employment+Wy_Employment,Total_Wages=US_Wages+Wy_Wages) %>% group_by(Period) %>% mutate(Rank=rank(-rank(Total_Employment))) %>% arrange(Period,Rank) %>% ungroup %>% mutate(Rank=ifelse(Rank>15,16,Rank))


RANK  <- EMP_SUMMARY %>% select(Occupation,Period,Rank)    %>% ungroup 
PEAK_WY <- RES %>% mutate(Period=ifelse(Year>=2030,"Operations","Construction")) %>% left_join(RANK) %>% mutate(Occupation=ifelse(Rank==16,"Other",Occupation)) %>% group_by(Occupation,Year,Period,Rank) %>% summarize(Peak_Employment=max(Employment,na.rm=TRUE)) %>% group_by(Period,Occupation,Rank)  %>% summarize(Wy_Peak_Employment=max(Peak_Employment,na.rm=TRUE)) %>% ungroup %>% select(Occupation,Period,Wy_Peak_Employment)

PEAK_US <- 	RES_US %>% select(Occupation,Year,US_Employment=Employment) %>% left_join(RES %>% mutate(Wy_Employment=Employment)) %>% mutate(Employment=US_Employment-Wy_Employment)%>% select(-US_Employment,-Wages) %>% mutate(Period=ifelse(Year>=2030,"Operations","Construction")) %>% left_join(RANK) %>% mutate(Occupation=ifelse(Rank==16,"Other",Occupation)) %>% group_by(Occupation,Year,Period,Rank) %>% summarize(Peak_Employment=max(Employment,na.rm=TRUE)) %>% group_by(Period,Occupation,Rank)  %>% summarize(US_Peak_Employment=max(Peak_Employment,na.rm=TRUE)) %>% ungroup %>% select(Occupation,Period,US_Peak_Employment)
rbind(RES_US,RES) %>% group_by(Year) %>% summarize(Employment=max(Employment)) %>% plot
PEAK_ALL <-	rbind(RES_US,RES) %>% group_by(Occupation,Year) %>% summarize(Employment=max(Employment)) %>% mutate(Period=ifelse(Year>=2030,"Operations","Construction")) %>% left_join(RANK) %>% mutate(Occupation=ifelse(Rank==16,"Other",Occupation)) %>% group_by(Occupation,Year,Period,Rank) %>% summarize(Peak_Employment=max(Employment,na.rm=TRUE)) %>% group_by(Occupation,Period) %>% summarize(Total_Peak_Employment=max(Peak_Employment)) %>% ungroup

#PEAK_ALL <- 	rbind(US ,WY) %>% group_by(Occupation,Year,Period) %>% summarize(Employment=sum(Employment,na.rm=TRUE)) %>% left_join(RANK) %>% mutate(Occupation=ifelse(Rank==16,"Other",Occupation)) %>% group_by(Occupation,Period) %>% summarize(Total_Peak_Employment=max(Employment,na.rm=TRUE)) %>% ungroup

PEAK <- PEAK_ALL  %>% left_join(PEAK_WY) %>% left_join(PEAK_US)
PEAK %>% arrange(desc(Total_Peak_Employment))
PEAK 


EMP_SUMMARY <- EMP_SUMMARY %>% select(-Peak_Wy_Employment,-Peak_US_Employment) %>% mutate(Occupation=ifelse(Rank>15,"Other",Occupation),Rank=ifelse(Rank>15,16,Rank)) %>% group_by(Occupation,Rank,Period) %>% summarize(Wy_Employment=sum(Wy_Employment,na.rm=TRUE),Wy_Wages=sum(Wy_Wages,na.rm=TRUE),US_Employment=sum(US_Employment,na.rm=TRUE),US_Wages=sum(US_Wages,na.rm=TRUE),Total_Employment=sum(Total_Employment,na.rm=TRUE),Total_Wages=sum(Total_Wages,na.rm=TRUE)) %>% ungroup
EMP_SUMMARY <- EMP_SUMMARY %>% left_join(PEAK) %>% arrange(Period,Rank) %>% print(n=100)

CONSTRUCTION_EMP_SUMMARY <- EMP_SUMMARY %>% select(Occupation,Period,Wy_Employment,Wy_Peak_Employment,Wy_Wages,US_Employment,US_Peak_Employment,US_Wages,Total_Employment,Total_Peak_Employment,Total_Wages) %>% filter(Period=='Construction') %>% select(-Period)

OP_EMP_SUMMARY <- EMP_SUMMARY %>% select(Occupation,Period,Wy_Employment,Wy_Peak_Employment,Wy_Wages,US_Employment,US_Peak_Employment,US_Wages,Total_Employment,Total_Peak_Employment,Total_Wages) %>% filter(Period=='Operations') %>% select(-Period)



############Save results
HEADER <- "Data Center_"
write.csv(EMP_MAIN_SUMMARY ,paste0("./Results/",HEADER,"Employment_Summary.csv"),row.names=FALSE)
write.csv(WY_NPV,paste0("./Results/",HEADER,"Wyoming_NPV_Values.csv"),row.names=FALSE)
write.csv(US_NPV,paste0("./Results/",HEADER,"US_NPV_Values.csv"),row.names=FALSE)
write.csv(NPV,paste0("./Results/",HEADER,"NPV_Values.csv"),row.names=FALSE)
write.csv(ALL_TAX,paste0("./Results/",HEADER,"Tax_Values.csv"),row.names=FALSE)
write.csv(CONSTRUCTION_EMP_SUMMARY ,paste0("./Results/",HEADER,"Construction_Occupation.csv"),row.names=FALSE)
write.csv(OP_EMP_SUMMARY,paste0("./Results/",HEADER,"Operation_Occupation.csv"),row.names=FALSE)


ggsave( filename = paste0("./Results/",HEADER,"Dollar Values of Model.png"),  plot = IMPLAN_REGION_DOLLAR_PLOT,  width = 8.5,  height = 11*(3./4),  units = "in", dpi = 300)
ggsave( filename = paste0("./Results/",HEADER,"Tax Model.png"),  plot = TAX_PLOT,  width = 8.5,  height = 11*(3./4),  units = "in", dpi = 300)
ggsave( filename = paste0("./Results/",HEADER,"Employment.png"),  plot = IMPLAN_EMPLOYMENT_PLOT ,  width = 8.5,  height = 11*(2/4),  units = "in", dpi = 600)