--- disqus: ahb0222 GA : G-VF9ZT413CG --- # 以R語言爬取監測站歷史資料並以ggplot2繪製風玫瑰圖(風花圖,Wind Rose)_大寮測站為例 > [color=#40f1ef][name=LHB阿好伯, 2021/02/15][:earth_africa:](https://www.facebook.com/LHB0222/) ###### tags: `R` `ggplot2` [TOC] ![](https://i.imgur.com/WRfddTo.jpg) :::danger 2024/01/20更新 非常可惜的是以前的網站更新了QQ 沒辦法在像之前一樣更改網址爬取資料 目前解決方式只好先分成兩步驟 先使用Python的Selenium下載CSV檔 再用ggplot2畫圖 相關程式碼我就分享出來 拋磚引玉看看有沒有大神能夠將整個程式碼在優化一下 或是有更好用的爬蟲套件能分享一下 不然光Selenium的安裝跟使用就耗費我好多時間XD Chrome的版本太新就是一個大問題chromedriver支援還沒跟上QQ 要想像之前架設一個Shiny能自動爬取跟畫圖看來是有困難QQ [使用Python selenium 爬取codis氣候資料服務系統資料並使用ggplot2繪製風花圖](/LW3Ys5qwQ4a4j6zrOdc8HA) ::: 今天來分享一個環工系比較可能用到的功能~~ 或是還有哪些科系會用到我就不清楚了 下面網站可以查詢台灣的監測資料 我們可以在上面查到許多台灣的降雨、風向、風速、日照等資料 ![](https://i.imgur.com/cpWgicg.jpg) 這次以公司附近的大寮測站為例 之前也有一篇相關文章大家可以看看[使用R(RStudio Cloud)語言擷取觀測資料查詢系統氣象資料](/cAVB2a4XRsC4SvAFhmygNQ) ![](https://i.imgur.com/zLH9aIJ.jpg) ![](https://i.imgur.com/1W1mcQZ.png) ```r= packages <- c("jsonlite", "rvest", "magrittr", "lubridate","ggplot2") installed_packages <- packages %in% rownames(installed.packages()) if (any(installed_packages == FALSE)) { install.packages(packages[!installed_packages]) } # 載入所需套件 invisible(lapply(packages, library, character.only = TRUE)) ``` :::danger ~~複製網址後刪除最後的時間資料 https://e-service.cwb.gov.tw/HistoryDataQuery/D..... &datepicker= 2021-05-01 貼到程式碼中 url_start <- "https://e-service.cwb.gov.tw/HistoryDataQuery/.....&datepicker=" 修改需要的時間範圍執行後續程式碼即可 start_date <- ymd("2021-01-01") #起始時間 end_date <- ymd("2021-03-31") #終止時間~~ 網站又...又...更新了 https://opendata.cwa.gov.tw/dataset/observation/O-A0001-001 https://codis.cwa.gov.tw/ 有空來調整看看 ::: ```r= #貼上網址刪除最後面的日期 #https://e-service.cwb.gov.tw/....&datepicker=2021-05-01 #https://e-service.cwb.gov.tw/....&datepicker= url_start <- "https://e-service.cwb.gov.tw/HistoryDataQuery/DayDataController.do?command=viewMain&station=C0V730&stname=%25E5%25A4%25A7%25E5%25AF%25AE&datepicker=" start_date <- ymd("2021-01-01") #起始時間 end_date <- ymd("2021-03-31") #終止時間 { data1 <- data.frame() data2 <- data.frame() for(d in c(0:(end_date - start_date))){ url <- paste0(url_start,as.character(start_date + d),"#") #時間資料 time_H <- url %>% read_html() %>% html_nodes(xpath='//*[(@id = "MyTable")]//td[(((count(preceding-sibling::*) + 1) = 1) and parent::*)]') %>% html_text(trim = T) %>% as.numeric() #風速資料 Windspeed <- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 7) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() #風向 windDirection<- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 8) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() data1 <- cbind(as.character(start_date + d), time_H, Windspeed,windDirection) #合併資料 data2 <- rbind(data2, data1) } data3 <- na.omit(data2) #刪除包含NA值的列 data3$windDirection_N <- cut(as.numeric(data3$windDirection), breaks = c(0, 11.26, 33.76, 56.26, 78.76, 101.26, 127.76, 146.26, 168.76, 191.26, 213.76, 236.26, 258.76, 281.26, 303.76, 326.26, 348.75, 360), labels = c("N","NNE","NE","ENE","E","ESE","ES","SSE","S","SWS","SW","WSW","W","WNW","NW","NNW","N"), # labels = c("北","北東北","東北","東東北","東","東東南","東南","南東南","南","南西南","西南","西西南","西","西西北","西北","西北西","北"), include.lowest = TRUE) data3$Windspeed_N <- cut(as.numeric(data3$Windspeed), breaks =c(0, 2.1, 4.1, 6.1, Inf), labels = c("0~2", "2.1~4", "4.1~6", ">6.1"), right = F) data3$Windspeed <- as.numeric(data3$Windspeed) ggplot(data = data3, aes(x = windDirection_N, fill = Windspeed_N))+ geom_bar( position = position_stack(reverse = TRUE)) + theme_bw() + scale_fill_discrete(guide = guide_legend(reverse=TRUE), name = "Wind Speed (m/s)")+ coord_polar(start = -0.3) + xlab("") } ``` ![](https://i.imgur.com/4h3OQhE.png) # Shiny APP ```r= # # This is a Shiny web application. You can run the application by clicking # the 'Run App' button above. # # Find out more about building applications with Shiny here: # # http://shiny.rstudio.com/ # packages <- c("jsonlite", "rvest", "magrittr", "lubridate","ggplot2","shiny") #install.packages() library("shiny") library("jsonlite") library("rvest") library("magrittr") library("lubridate") library("ggplot2") library("shiny") WindData <- function(station_number,start_date,end_date){ url_start <- paste0("https://e-service.cwb.gov.tw/HistoryDataQuery/DayDataController.do?command=viewMain&station=",station_number,"&stname=%25E5%25A4%25A7%25E5%25AF%25AE&datepicker=") start_date <- ymd(start_date) #起始時間 end_date <- ymd(end_date) #終止時間 data1 <- data.frame() data2 <- data.frame() for(d in c(0:(end_date - start_date))){ url <- paste0(url_start,as.character(start_date + d),"#") #時間資料 time_H <- url %>% read_html() %>% html_nodes(xpath='//*[(@id = "MyTable")]//td[(((count(preceding-sibling::*) + 1) = 1) and parent::*)]') %>% html_text(trim = T) %>% as.numeric() #風速資料 Windspeed <- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 7) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() #風向 windDirection<- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 8) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() data1 <- cbind(as.character(start_date + d), time_H, Windspeed,windDirection) #合併資料 data2 <- rbind(data2, data1) } windD <- factor(c("N","NNE","NE","ENE","E","ESE","ES","SSE","S","SWS","SW","WSW","W","WNW","NW","NNW","N")) data3 <- na.omit(data2) #刪除包含NA值的列 data3$windDirection_N <- cut(as.numeric(data3$windDirection), breaks = c(0, 11.26, 33.76, 56.26, 78.76, 101.26, 127.76, 146.26, 168.76, 191.26, 213.76, 236.26, 258.76, 281.26, 303.76, 326.26, 348.75, 360), windD, include.lowest = TRUE) data3$Windspeed_N <- cut(as.numeric(data3$Windspeed), breaks =c(0, 2.1, 4.1, 6.1, Inf), labels = c("0~2", "2.1~4", "4.1~6", ">6.1"), right = F) data3$Windspeed <- as.numeric(data3$Windspeed) p1 <- ggplot(data = data3, aes(x = windDirection_N, fill = Windspeed_N))+ geom_bar( position = position_stack(reverse = TRUE)) + theme_bw() + scale_x_discrete(drop = FALSE, labels = windD) + scale_fill_discrete(guide = guide_legend(reverse=TRUE), name = "Wind Speed (m/s)")+ coord_polar(start = -0.2) + xlab("") + ggtitle(paste(station_number,"(",start_date,"~",end_date,")")) p1 } WindData2 <- function(station_number,start_date,end_date){ url_start <- paste0("https://e-service.cwb.gov.tw/HistoryDataQuery/MonthDataController.do?command=viewMain&station=",station_number,"&stname=%25E5%25A4%25A7%25E5%25AF%25AE&datepicker=") start_date <- floor_date(ymd(start_date),unit = "month") #起始時間 end_date <- floor_date(ymd(end_date),unit = "month") #終止時間 int <- interval(start_date, end_date) m1 <- time_length(int, "month") data1 <- data.frame() data2 <- data.frame() for(d in c(0:m1)){ m2 <- start_date + months(d) url <- paste0(url_start,as.character(paste0(year(m2),"-",ifelse(month(m2)>=10, month(m2),paste0(0,month(m2))))),"#") #時間資料 time_d <- url %>% read_html() %>% html_nodes(xpath='//*[(@id = "MyTable")]//td[(((count(preceding-sibling::*) + 1) = 1) and parent::*)]') %>% html_text(trim = T) %>% as.numeric() #風速資料 Windspeed <- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 17) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() #風向 windDirection<- url %>% read_html() %>% html_nodes(xpath='//td[(((count(preceding-sibling::*) + 1) = 18) and parent::*)]') %>% html_text(trim = T) %>% gsub("X","-9999", .) %>% as.numeric() data1 <- cbind(as.character(start_date + months(d)), time_d, Windspeed,windDirection) #合併資料 data2 <- rbind(data2, data1) } data2 windD <- factor(c("N","NNE","NE","ENE","E","ESE","ES","SSE","S","SWS","SW","WSW","W","WNW","NW","NNW","N")) data3 <- na.omit(data2) #刪除包含NA值的列 data3$windDirection_N <- cut(as.numeric(data3$windDirection), breaks = c(0, 11.26, 33.76, 56.26, 78.76, 101.26, 127.76, 146.26, 168.76, 191.26, 213.76, 236.26, 258.76, 281.26, 303.76, 326.26, 348.75, 360), windD, include.lowest = TRUE) data3$Windspeed_N <- cut(as.numeric(data3$Windspeed), breaks =c(0, 2.1, 4.1, 6.1, Inf), labels = c("0~2", "2.1~4", "4.1~6", ">6.1"), right = F) data3$Windspeed <- as.numeric(data3$Windspeed) p1 <- ggplot(data = data3, aes(x = windDirection_N, fill = Windspeed_N))+ geom_bar( position = position_stack(reverse = TRUE)) + theme_bw() + scale_x_discrete(drop = FALSE, labels = windD) + scale_fill_discrete(guide = guide_legend(reverse=TRUE), name = "Wind Speed (m/s)")+ coord_polar(start = -0.2) + xlab("") + ggtitle(paste(station_number,"(",paste0(year(start_date),"-",month(start_date)), "~",paste0(year(end_date),"-",month(end_date)),")")) p1 } # Define UI for application that draws a histogram ui <- fluidPage( titlePanel("以R語言爬取監測站歷史資料並以ggplot2繪製風玫瑰圖(風花圖,Wind Rose)"), # Default value is the date in client's time zone # Application title # Sidebar with a slider input for number of bins sidebarLayout(textInput ( "station_number" , " 監測站編號_預設大寮測站", value = "C0V730" ), sidebarPanel(dateRangeInput('dateRange', label = 'Date range input: yyyy-mm-dd', start = Sys.Date() -1, end = Sys.Date() -1), submitButton(text = "更新", icon = NULL, width = NULL), downloadButton("downloadData", "Download Plot") )), helpText("測站編號查詢","https://e-service.cwb.gov.tw/HistoryDataQuery/"), radioButtons( "mode", "日/月資料",choices = c("日","月")), # Show a plot of the generated distribution mainPanel( plotOutput("WindRosePlote", width = "600px", height = "400px") ) ) # Define server logic required to draw a histogram server <- function(input, output) { output$WindRosePlote <- renderPlot( if(input$mode == "日"){ WindData(input$station_number, input$dateRange[1], input$dateRange[2]) } else{ WindData2(input$station_number, input$dateRange[1], input$dateRange[2]) } ) output$downloadData <- downloadHandler( filename = function(){paste("input$downloadData",'.png',sep='')}, content = function(file){ ggsave(file,plot=WindData(input$station_number, input$dateRange[1], input$dateRange[2])) } ) } # Run the application shinyApp(ui = ui, server = server) ``` ## [測試網站](https://gtgrthrst.shinyapps.io/Wind_Rose/) :::success * 修改風花圖有些角度沒有數據時不顯示的問題 * 解決方式在原本ggplot2的程式碼加上scale_x_discrete(drop = FALSE, labels = windD) * 增加"更新"按鍵減少原本更換時間就自動跑程式的問題 * 在shiny ui()中增加 submitButton(text = "更新", icon = NULL, width = NULL) 即可新增一個更新按鈕 * 增加日/月資料的選擇,若是需要一整季的資料可以選擇使用"月"資料減少等待時間 * 增加選擇按鈕radioButtons( "mode", "日/月資料",choices = c("日","月")) 寫兩個函數一個為使用日報表一個使用月報表資料 在server()中再判斷選擇的是什麼而執行哪個函數 ::: ![](https://i.imgur.com/gnRUNBV.png) # RStudio Cloud 專案 https://rstudio.cloud/spaces/104923/project/2596118 # [Shinydashboard(儀錶板)手動添加數據繪製風花圖](/Kd2uUuldTrmfafDFt_GfCA) ![](https://i.imgur.com/KFdMCHb.png) 🌟全文可以至下方連結觀看或是補充 全文分享至 https://www.facebook.com/LHB0222/ https://www.instagram.com/ahb0222/ 有疑問想討論的都歡迎於下方留言 喜歡的幫我分享給所有的朋友 \o/ 有所錯誤歡迎指教 # [:page_with_curl: 全部文章列表](https://hackmd.io/@LHB-0222/AllWritings) ![](https://i.imgur.com/47HlvGH.png)