星空研究Memo

ここは某天文屋の外部記憶装置である。

AIP4Win (.txt) -> Light curve

2011-03-15 19:24:27 | R言語

AIPのテキストから、一気に光度曲線かくコードを書いてみた。
主に連続測光時のクイックルック用に使う予定。
横軸はUTで、V-Cのエラーバーと、凡例も出るようにした。
凡例の位置を自動で決められるようにしたいが、プロットの具合によって凡例の位置は
ある程度微調整ができたほうが良いかな?!... まぁ、とりあえずできたコードを以下にメモ。

-----

txt_name <- "ER_UMa_20110304.txt"     # 入力名
png_name <- "ERUMA_20110304.png"      # 出力名
t <- "ER UMa (2011/03/04, by OUS)"    # グラフのタイトル
kc <- "K-C -1.1"                      # K-Cのラベル
ofs <- -1.1                           # K-Cのオフセット量
x <- 4.64                             # 凡例の位置x
y <- 0.57                             # 凡例の位置y

d <- read.table(txt_name, skip=38)    # AIP Ver.2のtxt読み込み
  attach(d)

 date <- V2
 time <- V3

 V_C <- V7
 V_C_err <- V8

 K_C <- V9 + ofs                      
 K_C_err <- V10

    D <- substr(date, 9, 10)          
     d  <- as.integer(D)
    H  <- substr(time, 1, 2)          
     h  <- as.integer(H)
    M  <- substr(time, 4, 5)         
     m  <- as.integer(M)
    S  <- substr(time, 7, 8)          
     s  <- as.integer(S)

 UT <- (((((s/60)+m)/60)+h)/24) + d

ylim = range(c(
  V_C + V_C_err,
  V_C - V_C_err,
  K_C + K_C_err,
  K_C - K_C_err))

VC <- data.frame(UT, V_C)
KC <- data.frame(UT, K_C)

#-----プロットと出力
png(png_name, pointsize = 15,
    width = 600, height = 540)

plot(KC, type="p", col="blue", cex=0.5,
     pch=3, ylim = rev(ylim), xlab = "UT", ylab = "relative magnitude")

par(xpd=TRUE)

legend(x, y,
       c("V-C", kc), col = c("red", "blue"),
       pch = c(16, 3), cex=0.8, box.lty=0)

arrows(UT, V_C + V_C_err,
       UT, V_C - V_C_err,
       col = 8,
       lwd = 1, lty = 1, length = 0)

points(VC, col="red", pch=16)
title(main=t, cex.main = 1)

par(new=T)
dev.off()

最新の画像もっと見る

post a comment