Rで自作関数:plot_shadow_path

関数名: plot_shadow_path

垂直な棒が地面に落とす影の先端の軌跡を描画する関数です。

上が北・右が東となる平面(X-Y座標)に、棒の位置を原点 (0, 0) として影の動きをプロットします。

本ポストはこちらの続きです。下記コード中の関数 get_sun_position もこちらで確認してください。

Rで自作関数:get_sun_position
関数名: get_sun_position緯度、経度、日時を入力として受け取り、度(Degrees)単位で高度と方位(北=0°, 東=90°, 南=180°, 西=270°)を返す関数です。参照資料 : #' 太陽の高度と方位を計算する関数...
Rで自作関数:plot_sun_path
関数名: plot_sun_path指定した日の「日の出」から「日の入り」までの太陽の軌跡(方位と高度)を計算し、グラフとして描画する関数です。本ポストはこちらの続きです。下記コード中の関数 get_sun_position もこちらで確認...
# パッケージの読み込み
library(ggplot2)
library(ggrepel)

#' 棒の日影図(Shadow Path)をプロットする関数
#'
#' @param date 日付(文字列 "YYYY-MM-DD" または Date型)
#' @param lat 緯度(度、10進数)
#' @param lon 経度(度、10進数)
#' @param pole_height 垂直な棒の高さ
#' @param tz タイムゾーン(デフォルトは "Asia/Tokyo")
#' @return ggplotオブジェクト
plot_shadow_path <- function(date, lat, lon, pole_height, tz = "Asia/Tokyo") {
  
  target_date <- as.Date(date)
  deg2rad <- function(deg) deg * pi / 180
  
  # 1. 指定日の00:00から23:59まで「1分間隔」の時刻ベクトルを作成
  start_time <- as.POSIXct(sprintf("%s 00:00:00", target_date), tz = tz)
  end_time   <- as.POSIXct(sprintf("%s 23:59:00", target_date), tz = tz)
  time_seq_1min <- seq(start_time, end_time, by = "1 min")
  
  # 2. 太陽位置を一括計算
  sun_pos <- get_sun_position(lat = lat, lon = lon, datetime = time_seq_1min)
  
  df_all <- data.frame(
    datetime  = time_seq_1min,
    elevation = sun_pos$elevation,
    azimuth   = sun_pos$azimuth
  )
  
  # 3. 影が存在する条件(太陽高度が 0度より大きい)で抽出
  df_daylight <- subset(df_all, elevation > 0)
  if (nrow(df_daylight) == 0) stop("影が形成されるほど太陽が昇りません。")
  
  # 4. 影の長さとX-Y座標の計算
  df_daylight$shadow_len <- pole_height / tan(deg2rad(df_daylight$elevation))
  df_daylight$x <- -df_daylight$shadow_len * sin(deg2rad(df_daylight$azimuth))
  df_daylight$y <- -df_daylight$shadow_len * cos(deg2rad(df_daylight$azimuth))
  
  # 5. プロット用にデータを間引く(5分間隔)
  is_5min_interval <- as.numeric(format(df_daylight$datetime, "%M")) %% 5 == 0
  df_plot <- df_daylight[is_5min_interval, ]
  
  # 時刻ラベル用(毎正時 00分)のデータ
  is_on_the_hour <- format(df_plot$datetime, "%M") == "00"
  df_labels <- df_plot[is_on_the_hour, ]
  
  # --- グラフの表示範囲を固定 ---
  plot_limit <- ceiling(max(abs(df_labels$x),abs(df_labels$y)))

  # 6. ggplot2 で日影図を描画
  p <- ggplot(df_plot, aes(x = x, y = y)) +
    
    # 基準線(東西・南北)
    geom_hline(yintercept = 0, color = "gray70", linewidth = 0.5, linetype = "dashed") +
    geom_vline(xintercept = 0, color = "gray70", linewidth = 0.5, linetype = "dashed") +
    
    # 影の軌跡
    geom_path(aes(color = as.numeric(datetime)), linewidth = 1.5, show.legend = FALSE) +
    scale_color_viridis_c(option = "C") +
    
    # 原点(垂直な棒の位置)
    geom_point(aes(x = 0, y = 0), size = 4, color = "black", shape = 15) +
    
    # 時刻のポイント
    geom_point(data = df_labels, size = 3, color = "red") +
    
    # --- ggrepelでラベルの重なりを回避 ---
    geom_text_repel(
      data = df_labels, 
      aes(label = format(datetime, "%H:%M")), 
      size = 3.5, 
      color = "black",
      box.padding = 0.6,       # ラベル同士の余白
      point.padding = 0.3,     # ポイントとラベルの余白
      segment.color = "gray50",# 引き出し線の色
      max.overlaps = Inf       # 全てのラベルを強制表示
    ) +
    
    # --- coord_fixedに limits を設定してズーム ---
    coord_fixed(
      xlim = c(-plot_limit, plot_limit), 
      ylim = c(-plot_limit, plot_limit)
    ) +
    
    # タイトルと軸ラベル
    labs(
      title = sprintf("日影曲線(Shadow Path)- %s", target_date),
      subtitle = sprintf("緯度: %.4f, 経度: %.4f / 棒の高さ: %.1f\n■: 棒の位置", 
                         lat, lon, pole_height),
      x = "西 (-x)   ←   東西方向の影の長さ   →   東 (+x)",
      y = "南 (-y)   ←   南北方向の影の長さ   →   北 (+y)"
    ) +
    
    theme_bw(base_size = 12) +
    theme(
      plot.title = element_text(face = "bold", size = 16),
      panel.grid.minor = element_blank()
    )
  return(p)
}

実行例

# 東京の緯度経度
my_lat <- 35.689
my_lon <- 139.691

# 棒の高さ
my_pole_height <- 10.0 

# 夏至(2026年6月21日)のプロット
plot_summer_shadow <- plot_shadow_path(date = "2026-06-21",lat =  my_lat,lon =  my_lon,pole_height =  my_pole_height)
print(plot_summer_shadow)

# 冬至(2026年12月22日)のプロット
plot_winter_shadow <- plot_shadow_path("2026-12-22", my_lat, my_lon, my_pole_height)
print(plot_winter_shadow)
Figure 1
Figure 2

以上です。