library(rvest)
library(lubridate)
library(dplyr)
library(ggplot2)
library(gridExtra)

setwd("~/RProjects/jakosc_powietrza")

# czy dane ze strony czy z przygotowanego pliku?
online <- FALSE
# ile ostatnich dni wczytać?
ile_dni <- 1095 # 3 lata

if(online == TRUE) {
   # wczytaj dane ze strony
   # dzisiaj
   dzis <- Sys.Date()
   dzis <- make_date(year(dzis), month(dzis), day(dzis))
   
   # dane bazowe
   sojp_base_url <- "http://sojp.wios.warszawa.pl/?page=raport-dobowy&data="
   dane_dobowe <- data.frame()

   # pobieraj kolejne dni
   for(i in 1:ile_dni) {
   
      # progress bar :)
      cat(paste("\rPobieram dane:", i, "/", ile_dni, "\r"))
   
      # kolejne dni wstecz, począwszy od wczoraj
      dzien <- dzis-ddays(i)
   
      # przygotowanie urla dla konkretnej daty
      strona <- read_html(paste0(sojp_base_url, format(dzien, "%d-%m-%Y")))
   
      # pobranie danych do tabeli
      df <- strona %>% html_node("table") %>% html_table(fill=TRUE)
   
      # w pierwszym wierszu danych tabeli są jednostki - usuwamy
      df <- df[-1,]
   
      # konwersja na liczby
      for(j in 2:13) {
         df[,j] <- as.numeric(df[,j])
      }
   
      # dodajemy datę pomiaru
      df <- cbind(df, data=dzien)
   
      # porządek w nazwach kolumn - jednostki µg/m3 + %LV (% dopuszczalnego w Polsce poziomu)
      colnames(df) <- c("Stacja", "PM10", "PM10_pLV", "NO2", "NO2_pLV",
                        "CO", "CO_pLV", "O3", "O3_pLV", "SO2_S1", "SO2_S1_pLV",
                        "SO2_S24", "SO2_S24_pLV", "Data")
   
      # scalenie ze wszystkimi
      dane_dobowe <- rbind(dane_dobowe, df)
   }
   
   # trafiają się wiersze puste - z pustą nazwą stacji -> usuwamy
   dane_dobowe <- dane_dobowe[nchar(dane_dobowe$Stacja)!=0,]

   # zmiana wartości liczbowych na opisowe
   # skale zgodne z opisem na http://sojp.wios.warszawa.pl/
   dane_dobowe$PM10_jakosc <- cut(dane_dobowe$PM10,
                                  breaks = c(0,20,60,100,140,200,1e+6),
                                  labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   dane_dobowe$NO2_jakosc <- cut(dane_dobowe$NO2,
                                 breaks = c(0,40,100,150,200,400,1e+6),
                                 labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   dane_dobowe$CO_jakosc <- cut(dane_dobowe$CO,
                                breaks = c(0,2,6,10,15,20,1e+6),
                                labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   dane_dobowe$O3_jakosc <- cut(dane_dobowe$O3,
                                breaks = c(0,30,70,120,160,240,1e+6),
                                labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   dane_dobowe$SO2_S1_jakosc <- cut(dane_dobowe$SO2_S1,
                                    breaks = c(0,50,100,200,350,500,1e+6),
                                    labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   dane_dobowe$SO2_S24_jakosc <- cut(dane_dobowe$SO2_S24,
                                     breaks = c(0,5,10,15,20,50,1e+6),
                                     labels = c("Bardzo dobry", "Dobry", "Umiarkowany", "Dostateczny", "Zły", "Bardzo zły"))
   
   # zapisanie danych do pliku - na wszelki wypadek,
   # żeby nie pobierać za każdym razem w przyszłości
   save(dane_dobowe, file="dane_dobowe.Rdata")
} else {
   # wczytaj dane z pliku
   load("dane_dobowe.Rdata")   
}


# lista wszystkich stacji
stacje <- unique(dane_dobowe$Stacja)
print(stacje)
lp Stacja
1 Belsk-IGFPAN
2 Granica-KPN
3 Guty Duże
4 Legionowo-Zegrzyńska
5 Otwock-Brzozowa
6 Piastów-Pułaskiego
7 Płock-Reja
8 Płock-Gimnazjum
9 Radom-Tochtermana
10 Siedlce-Konarskiego
11 Warszawa-Komunikacyjna
12 Warszawa-Marszałkowska
13 Warszawa-Podleśna
14 Warszawa-Targówek
15 Warszawa-Ursynów
16 Żyrardów-Roosevelta 
# wybrana stacja
stacja <- "Warszawa-Targówek"


wykres <- function(dane, wybrana_stacja, wskaznik_nazwa, wskaznik_jakosc, tytul) {
   # funkcja przygotowuje wykres - na podstawie wskazanego wskaźnika, z podanym tytułem
   
   # pomocnicze dane
   dane_df <- data.frame(Data=dane$Data,
                         Stacja=dane$Stacja,
                         wskaznik=dane[[wskaznik_nazwa]],
                         wskaznik_jakosc=dane[[wskaznik_jakosc]])

   # jeśli są dane - narysuj wykres
   # jeśli nie ma - zrób zaślepkę
   if(dane_df %>% filter(Stacja==wybrana_stacja, !is.na(wskaznik)) %>% nrow() != 0) {
      p <- dane_df %>%
         filter(Stacja==stacja, !is.na(wskaznik)) %>%
         ggplot() +
         geom_point(aes(Data, wskaznik, color=wskaznik_jakosc)) +
         theme_minimal() +
         labs(title=tytul, x="Data", y=wskaznik_nazwa, colour="Jakość") +
         ylim(0,max(dane_df$wskaznik, na.rm = TRUE)) +
         scale_color_brewer(palette="Set1", direction = -1)
   } else {
      p <- ggplot() +
         theme_minimal() +
         labs(title=paste0(tytul, "\n\nBRAK DANYCH ZE STACJI"))
   }
   
   return(p)
}

# przygotowanie wykresów wszystkich wskaźników
p1 <- wykres(dane_dobowe, stacja, "NO2_pLV", "NO2_jakosc", "NO2 - maksymalna średnia 1h")
p2 <- wykres(dane_dobowe, stacja, "PM10_pLV", "PM10_jakosc", "Pył zawieszony PM10 - średnia 24h")
p3 <- wykres(dane_dobowe, stacja, "CO_pLV", "CO_jakosc", "CO - maksymalna średnia 8")
p4 <- wykres(dane_dobowe, stacja, "O3_pLV", "O3_jakosc", "O3 - maksymalna średnia 8h")
p5 <- wykres(dane_dobowe, stacja, "SO2_S1_pLV", "SO2_S1_jakosc", "SO2 - maksymalna średnia 1h")
p6 <- wykres(dane_dobowe, stacja, "SO2_S24_pLV", "SO2_S24_jakosc", "SO2 - średnia 24h")
# wszystkie wykresy na raz - prezentacja
grid.arrange(p1, p2, p3, p4, p5, p6,
             top=paste("Procent dopuszczalnej normy dla wskaźników\nna stacji", stacja))