2018-08-25

VMAFで画質比較、パラメータ検討

DTV

h264_vaapiでハードウェアエンコードしてたんだが、素材に寄ってはブロックノイズが気になることがあり、もう一度、コーデックも含め、パラメータを検討してみることにした。目視では時間がかかるうえ、客観性も劣る。そこで、VMAFという指標を用いて比較することにした。

VMAFはffmpegのフィルタとして実装されており、エンコ前後の動画を食わせると100点満点で値を返してくる。手元のffmpegをコンパイルし直すのが面倒だったので、ニコニコで配ってるバイナリを引っ張ってきて使用した。

比較するコーデックはh264, h264_vaapi, hevc、それぞれ固定画質モードでエンコードし、ファクタは17から33までひとつおきとした。素材は4種類。洋ドラ、演芸、動物物、アニメ。こんな感じ。

ffmpeg -loglevel 32 -hide_banner -analyzeduration 30M -probesize 30M -i CUT.ts -map 0:0 -map 0:1 -vcodec $CODEC -flags +ildct+ilme -crf $CRF -acodec copy -f mpegts -mpegts_m2ts_mode 1 -vsync 0 $OUTFILE

vaapiは無いな、と。わざわざエンコードする意味がない。HWエンコはやはり配信用だ。
hevcとh264、いい勝負なんだが、素材によって得手不得手が有るようだ。
アニメはh264がいい仕事している。
動きの少ない演芸番組、洋ドラ、低ビットレートではhevcが勝る。
VMAF90点を目標とするなら、hevcでもh264でもqp=27としてエンコードするのがバランスが良いだろう。
エンコード時間を考えると、h264だな。
改めてh264の優秀さに脱帽。

#!/usr/bin/R

vmaf <- NULL
vmaf <- data.frame (
  sample=c("endeavour", "endeavour", "endeavour", "endeavour", "endeavour",
           "endeavour", "endeavour", "endeavour", "endeavour",
           "endeavour", "endeavour", "endeavour", "endeavour", "endeavour",
           "endeavour", "endeavour", "endeavour", "endeavour",
           "endeavour", "endeavour", "endeavour", "endeavour", "endeavour",
           "endeavour", "endeavour", "endeavour", "endeavour",

           "engei", "engei", "engei", "engei", "engei",
           "engei", "engei", "engei", "engei",
           "engei", "engei", "engei", "engei", "engei",
           "engei", "engei", "engei", "engei",
           "engei", "engei", "engei", "engei", "engei",
           "engei", "engei", "engei", "engei",

           "neko", "neko", "neko", "neko", "neko", 
           "neko", "neko", "neko", "neko", 
           "neko", "neko", "neko", "neko", "neko", 
           "neko", "neko", "neko", "neko", 
           "neko", "neko", "neko", "neko", "neko", 
           "neko", "neko", "neko", "neko", 

           "baka", "baka", "baka", "baka", "baka", 
           "baka", "baka", "baka", "baka", 
           "baka", "baka", "baka", "baka", "baka", 
           "baka", "baka", "baka", "baka", 
           "baka", "baka", "baka", "baka", "baka", 
           "baka", "baka", "baka", "baka"), 

  codec=c("h264", "h264", "h264", "h264", "h264",
          "h264", "h264", "h264", "h264",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "hevc", "hevc", "hevc", "hevc", "hevc",
          "hevc", "hevc", "hevc", "hevc",

          "h264", "h264", "h264", "h264", "h264",
          "h264", "h264", "h264", "h264",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "hevc", "hevc", "hevc", "hevc", "hevc",
          "hevc", "hevc", "hevc", "hevc",

          "h264", "h264", "h264", "h264", "h264",
          "h264", "h264", "h264", "h264",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "hevc", "hevc", "hevc", "hevc", "hevc",
          "hevc", "hevc", "hevc", "hevc",

          "h264", "h264", "h264", "h264", "h264",
          "h264", "h264", "h264", "h264",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "h264_vaapi", "h264_vaapi", "h264_vaapi", "h264_vaapi",
          "hevc", "hevc", "hevc", "hevc", "hevc",
          "hevc", "hevc", "hevc", "hevc"),

  q=c(17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,

      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,

      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,

      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33,
      17, 19, 21, 23, 25, 27, 29, 31, 33),

  size=c(40233408, 28410816, 19067712, 13342848, 10179840,
         7934592, 6602880, 5611200, 4888128,
         58493376, 42685632, 34079040, 20838912, 17428800,
         12104064, 9887808, 8582208, 7333440,
         30403584, 19525248, 13114752, 9402624, 7192128,
         5808192, 4856448, 4212480, 3741504,

         151035264, 114034176, 82104768, 56709504, 39038208,
         26509056, 18812928, 13984704, 10893888,
         180936576, 129821568, 96960768, 60131904, 45524544,
         31225920, 22638912, 18226368, 14008320,
         84664320, 57035904, 37372224, 24697152, 16933632,
         12134592, 9122304, 7165824, 5833920,

         154437312, 122170176, 92806464, 67483776, 48272256,
         33981696, 25192320, 19374720, 15513600,
         195580608, 144600000, 111267264, 72640704, 56917248,
         40729728, 30924096, 25297920, 19798848,
         96467904, 71329152, 50634048, 34667328, 23605248,
         16603968, 12332160, 9608064, 7767552, 

         49433856, 38120256, 28001280, 20700864, 15813888,
         12406848, 10219392, 8692032, 7590144,
         89130048, 71908608, 61507008, 44960064, 39323328,
         30558528, 25833216, 23206080, 20117568,
         34869312, 26099328, 20038272, 15844416, 12858240,
         10714560, 9045696, 7749888, 6715392), 

  vmaf=c(97.414223, 97.033137, 96.482038, 95.757239, 94.927476,
         93.827896, 92.298622, 90.533203, 88.224027,
         97.030416, 96.358609, 95.820635, 94.789654, 93.996665,
         92.554915, 90.360646, 88.345397, 85.164701,
         96.219444, 95.691093, 95.029160, 94.137422, 92.967064,
         91.447516, 89.471573, 87.011807, 83.784730,

         98.938432, 98.642856, 98.124073, 97.340494, 96.344059,
         94.954042, 93.197262, 91.164133, 88.550417,
         98.414272, 97.745381, 97.169063, 95.961049, 94.971345,
         93.369077, 91.096610, 88.941571, 85.402841,
         96.982160, 96.112228, 95.127002, 94.008808, 92.675871,
         90.938491, 88.704656, 85.810726, 82.135652,

         98.956203, 98.674420, 98.200277, 97.477086, 96.427196,
         94.767823, 92.285201, 89.003958, 84.911819,
         98.546474, 97.897274, 97.294441, 95.725487, 94.356864,
         91.833858, 88.033021, 84.286361, 78.460347,
         97.610464, 96.862611, 95.852296, 94.467861, 92.648679,
         90.275109, 87.202078, 83.337630, 78.449135,

         97.403659, 97.187608, 96.859229, 96.446676, 96.002653,
         95.345043, 94.579849, 93.621566, 92.440076,
         97.229073, 96.942756, 96.600447, 96.138650, 95.621452,
         94.792517, 93.501437, 92.044247, 89.634072,
         96.659194, 96.318148, 95.893010, 95.349360, 94.573302,
         93.521138, 92.015785, 90.099163, 87.699080)
)
vmaf <- data.frame (vmaf, bitrate=vmaf$size / 1000 / 60 * 8)

png ("vmaf.png", width=960, height=540)

xr=c(500, 30000)
lt="x"

# 
  x <- subset (vmaf, sample=="endeavour"); cc=hsv (0, 1, 1);

  plot (x[which(x$codec=="h264"),]$bitrate, x[which(x$codec=="h264"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    xlab="Bitrate (kb/s)", ylab="VMAF Score", type="b", col=cc, pch=0)

  par(new=T)

  plot (x[which(x$codec=="h264_vaapi"),]$bitrate,
    x[which(x$codec=="h264_vaapi"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=1)

  par(new=T)

  plot (x[which(x$codec=="hevc"),]$bitrate, x[which(x$codec=="hevc"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=2)

  par(new=T)

#
  x <- subset (vmaf, sample=="engei"); cc=hsv (0.25, 1, 1)
 
  plot (x[which(x$codec=="h264"),]$bitrate, x[which(x$codec=="h264"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=0)

  par(new=T)

  plot (x[which(x$codec=="h264_vaapi"),]$bitrate,
    x[which(x$codec=="h264_vaapi"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=1)

  par(new=T)

  plot (x[which(x$codec=="hevc"),]$bitrate, x[which(x$codec=="hevc"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=2)

  par(new=T)

#
  x <- subset (vmaf, sample=="neko"); cc=hsv (0.5, 1, 1)

  plot (x[which(x$codec=="h264"),]$bitrate, x[which(x$codec=="h264"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=0)

  par(new=T)

  plot (x[which(x$codec=="h264_vaapi"),]$bitrate,
    x[which(x$codec=="h264_vaapi"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=1)

  par(new=T)

  plot (x[which(x$codec=="hevc"),]$bitrate, x[which(x$codec=="hevc"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=2)

  xx <- subset (x, codec=="h264_vaapi")
  text (xx$bitrate, xx$vmaf, xx$q, pos=4, col=cc, cex=0.75)
  # text (x[9,]$bitrate, x[9,]$vmaf, x[9,]$codec, pos=4, col=cc)
  par(new=T)

#
  x <- subset (vmaf, sample=="baka"); cc=hsv (0.75, 1, 1)

  plot (x[which(x$codec=="h264"),]$bitrate, x[which(x$codec=="h264"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=0)

  par(new=T)

  plot (x[which(x$codec=="h264_vaapi"),]$bitrate,
    x[which(x$codec=="h264_vaapi"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=1)

  par(new=T)

  plot (x[which(x$codec=="hevc"),]$bitrate, x[which(x$codec=="hevc"),]$vmaf,
    log=lt, xlim=xr, ylim=c(70, 100),
    axes=F, xlab="", ylab="", type="b", col=cc, pch=2)

  legend (10000, 90, c("h264", "h264_vaapi", "hevc"), pch=0:3)
  legend (10000, 84, c("endeavour", "engei", "neko", "baka"), col=c(hsv(c(0, 0.25, 0.5, 0.75), 1, 1)), lty=1)
  title ("VMAF Comparison")

dev.off()
#

2018-08-20

きょうの料理

最近の「きょうの料理ビギナーズ」、汚鮮されてね？お好み焼きと称してパジョンだったり、焼き肉のタレと称してプルコギだったり。酷えな。〇ね。〇〇人。藤野嘉子、おまえチョンだら。

2018-08-20

サマータイムだと

技術的にどうなのか、とか、人間工学的にどうなのか、とか、いろいろ言われているようだ。もち、ワシも反対だ。

それら反対者の意見に一言付け加えさせてもらいたい。言い出しっぺは単に馬鹿なのではなくて国民を馬鹿にしているんだ、と。そもそも、サマータイムって、「明日から一時間早く出勤してねー、関係公共機関もよろしくねー」ってお願いしても対応がまちまちで音頭が取れないから、「時計いじっちゃうし、おめーら黙って従え」ってことだと思うんだ。ゴリゴリのゴリ押しだわ。

そもそも政府に時間をいじる権利なんてあるのか。無いと思うけどなあ。

2018-08-15

パソコンが静かになりました

PCで録画している都合上、夜中に自動起動したりする。最近、音が気になってきてのでCPUファン、ケースファン共々交換した。現在使用しているケースは十年ほど前に購入したケース、CPUは数年前ぐらいかな。CPUファンは良しとしても、ケースファンの方は少なく見積もっても一万時間は動作してるよな。

購入したファンはこいつら、NOCTUAだ。静音・劇冷えだっ。ケース前面の80mm角ファンから吸気 → HDDを冷却 → CPUのヒートシンクへ → ケース後面の120mm角ファンで排気、という流れ。ヒートシンクが重くて怖いので、ケースは横置きにした。*1

[Noctua正規代理店]NH-D15S - サイレントタワーヒートシンク[NH-D15S]

出版社/メーカー: Noctua
メディア: Personal Computers
この商品を含むブログ (1件) を見る

Noctua - NF-S12A PWM対応ファン

出版社/メーカー: noctua
メディア: Personal Computers
この商品を含むブログ (1件) を見る

［Noctua正規販売代理店］NF-A8 PWM - 80mm Premium Quality Fan [NF-A8 PWM]

出版社/メーカー: Noctua
メディア: Personal Computers
この商品を含むブログを見る

すっげー静かになった。うるさい周波数で-15dBA、そうでもないところで-10dBAくらい。感覚的にはベランダのサッシ閉めたくらいか。アイドル時のCPU温度も5°C以上下がった。

ただ、1400Hzの音が耳障りといえば耳障りだ。恐らく5400rpmのHDDのケースへの共振だと思われる。スピンドルモータが12極で、5400÷12=450、で三次の共振で1350Hzってな感じでしょう。

しかも交換後のほうが鋭いピーク、多分ケースを掃除したとかネジの締め付け具合で共振の具合が変わったんだと思う。つまり、制振する手立てはありそうだということだ。なんとかせねば。

*1:というか、ヒートパイプの動作原理から言うと気化側が下、凝集側が上、なんだわ。ウィックの毛細管現象で重力方向を克服できるなどというが、あまり効果がなかったことを今になって思い出した。

2018-08-13

アンセキュア・マンション、D-Room

マンションの裏口、自転車置き場にこういうのがぶら下がってた。ま、あれだ、小箱付き南京錠ってやつだ。５桁くらいの数字で開閉し、なかにモノを入れて置けるっつうやつ。これが、マンションの敷地内とは言え、誰でも入れる場所にぶら下げてあるわけだ。民泊の運用で、事務所の位置を知られたくないなどの理由で使われることが多い。南京錠の中に、部屋の鍵とかカードキーとか入れといて宿泊客へ受け渡しをするわけだ。

気になったので管理会社に問い合わせたらこんな回答だった。「民泊ではありません。出入りの業者さんが利用するカードキーをそこに入れてある。業者さんが一旦大和リビングの事務所までカードキーを取りに来なくても良いようにそうしている。」んー、困ったぞ、こいつ馬鹿だぞ。カードキーはFelicaと同じ暗号強度があるっつうんで、ダイワハウスのCMでセキュア・マンションなどと謳ってるんだろよ。トリプルDESの相互認証を整数五桁までイッキにコンプロマイズする大胆さに驚き。指摘しても知らん顔の馬鹿さに二度驚き。*1 *2

こんなこともあった、これもカードキー絡みなんだが。大和リビングには錠と鍵の区別が出来ない奴が居るらしい。こういう掲示物は上長のチェックが入るのが普通だと思うが、チェックした人も区別の出来無い人なのか、そもそもチェックなんてしてないのか。どっちも嫌だね。その上、日本語も怪しい。ここに貼っとくから読んでみてね、こっちまで混乱してくるから。

*1:大和ハウスのHPから問い合わせをしたくリンクを辿るが、賃貸物件は全て大和リビングへの問い合わせに行き着く。しかたないので、大和ハウスの分譲物件の問い合わせメール窓口に、こうこうこういうことなんで、どうにかして下さい、とお願いしておいた。二日後、南京錠は撤去された。セキュア・マンション再び。しかし、問題は根が深くて、馬鹿が管理していてそれをチェックする人間が居ないっつうことじゃねえかと。市内に大和の賃貸物件が雨後の筍なんだが、管理に手が回らなくなった大和リビングはまた同様のポカをやると思う。

*2:その三ヶ月ほどしたら別の場所に南京錠がぶら下がってました。バカの相手はだるいので近々消費者センターにねじ込みます。

2018-06-15

Alien ^3

映画

BSでエイリアン一挙放送なんてやってるんでつい見てしまった。で、歳を取ると映画の見方が変わるもんだなと思った。しっとり古典の第一作、どんぱちバイオレンスの第二作、当時少年だった私はどちらも大いに楽しんだ。そして第三作目、なんか雰囲気は良くて抗いがたいものがあるが、噺がどーも良く分からん、と、好きになりたいけど好きになれねーな、みたいな。そんなもんでしたね。

そらそーだわ、中学生に理解できる映画ではなかった。ただそれだけのことでも無いことが、今回見直してみて判明した。そもそも訳を付けとる人が理解してない。そんな糞字幕で噺が理解できるわけ無い。今回のスターチャンネルでの放送、戸田奈津子を超える糞っぷり、見てみましょう。

会話が理解できない糞っぷり

字幕が酷くて見るのやめるレベル、つか止めた。で、今これ書いてる。ニュートの解剖、食堂での一幕の後、チャールズ・ダンス演じる医官とシガニー・ウィーバー演じるリプリーが身の上話をするシーン。

[Ripley] How did you get this wonderful assignment?
[Clemens] How do you like your new haircut?
[Ripley] It's okay.

「こんなクソみたいな仕事、志願したの？」と興味半分、ディスり半分で訊いてるわけよ。で、この小娘め、んなわけねーだろ、と質問に質問で返答する。「新しいヘアスタイルは気に入ったかい？」と。気に入るわけねえわな、剃られてんだもん。「しょうがないものね。」と返事しながら、医官の身の上を自分に重ねてんだよ。新しいヘアカットや靴、帽子などを職や立場の比喩とするのは欧米の会話では散見される。そうであることを知らずとも、この短い会話の後にやや長めの無言のシーンがあることで、先ほどの会話の反芻をやってる感は伝わってくるはずだ。この訳をつけた奴は英語を文字通り理解はできているが、会話を理解できていない。馬鹿だ、周囲にいても関わりたくないタイプ。

この後も突然「私に惹かれてるの？」と尋ねるリプリーに対して「どういう意味だ？」「そういう意味よ。」という会話が有る。これはさっきの職業志願の会話のリプライズ。やられっぱなしじゃイヤのリプリーが仕掛けたレトリックだ。それを「セックチュ」と訳しやがった馬鹿、出てこい。アタマ沸いとんのかい。

単語を知らない糞っぷり

誤訳も有るぜ。前後するが、身の上話の前、所長や施設についての会話中、語れば長いしややこしいので、会話の途中で一杯勧めている。

[Clemens] Tincture?

これを「２錠」と訳しやがった馬鹿、やっぱり出てこい。*1

ダブルミーニング、掛詞を訳に落とし込めない低能っぷり

剃髪して食堂に現れたリプリー、性衝動を隠しきれない囚人、リーダー的存在のディロンが信仰の話に重ねて、リプリーに警告しているシーン。

[Dillon] You see, we've got a good place to wait here.
[Dillon] And until now, no temptation.
[Ripley] What are you waiting for?

「な、分かるだろ？なんせ待ち時間はたっぷり有る。」「これまでには無かった誘惑もある。」と語るディロン。生長して信仰を遂げるには誘惑に打ち勝つ意思が大切だと説くと同時に、「お前レイプされちゃうぜ、気ぃつけや。」とも言ってる。これに「一体、何を待ってるの？」と応えるリプリー。信仰のコンテクストでは「キリストの再来的なものを待ってるんでしょ？あたしゃぁ信じませんが」という意味にも取れるし、レイプ警報のコンテクストでは「やれるもんならやってみなさいよ」という意味だ。ダブルミーニングにダブルミーニングで返すクレバーさ、そして鼻っ柱の強さ。このシーンだけでもリプリーのキャラクターがよく表現されている。…筈なのに、この盆暗な字幕のせいで台無し。

作品への敬意が感じられない

上で挙げたシーンはエイリアン登場前の導入部。どうせ皆死ぬんだからテキトーで良くね？的な和訳にゲンナリ。作品への敬意も、映画への愛も、和訳する能力も欠けている方々。ほんっとに辞めて欲しい。

苦情はこの辺までとして

しかし、良い映画だわ、海外でも国内でも評価は別れたと記憶しているが。役者、演技、台詞が絶妙に咬み合っていて良い。

チャールズ・ダンス(イングランド)、ポール・マッギャン(イングランド)、ブライアン・グローバー(イングランド)、ラルフ・ブラウン(イングランド)、ダニー・ウェッブ(イングランド)、ホルト・マッキャラニー(イングランド)、ピート・ポスルスウェイト(イングランド)と、イングランドの名優がズラリ。かの国の大抵の役者さんは俳優学校でガッチリ勉強してるからか、セリフ回しや演技がゴイスー。しかもギャラが安いので制作側もウッハー。

このイングランド勢が労働監獄という雰囲気を、先ず見た目で醸し出す。米人の役者とは明らかに顔つきが違う。歯並びも米人に比べれば悪い。で、喋りで醸し出す。米人にしてみれば馴染みのない、というか異国のアクセント。

まあ、評価が別れるのはしょうがない。待ちに待ったエイリアン新作がこんなだったらひっくり返るわな。

*1:後日吹き替え版で確認したところ、翻訳は石原千麻だった。向いてないから職業替えなさい。

2018-06-12

貧者の赤道儀（The Poor Man's Equatorial） - (2-h) フィールドソルバの精度

写真パノラマ

天体写真を撮って、フィールドソルバに掛けたとき、どのくらいの精度が得られるのだろうか。これまでのような用途であれば、1ピクセル以下だったら実用上問題なさそうでは有る。そもそもそういった精度が期待できるのか考えてみる。

FITSの座標系については、下から辿って、PAPER Iを読むとよく分かる。

FITS WCS Page

前回のパノラマの検証

前回のパノラマ、うまく繋がったように見えてたけど結構綱渡りだったんだろうな、と思い、検証してみることにした。

solve-fieldを実行すると結果ファイルがいろいろできるんだが、このなかに*.corrというのがある。correlationつうことだ。同定に使用した画像上の天体と天体カタログ上の天体、両方の座標が並べて表にしてくれている。field_{x,y}つうのが画像上の天体のピクセル座標、index_{x,y}がカタログ天体のピクセル座標だ。こいつらの距離が画像上でのズレになる。

大体画像一枚につき20から30くらいの天体が*.corrファイルにあるので、平均や偏差で評価できそうだ。とある画像で計算してみるとσ≒0.007だったので3σで考えても0.02ピクセルとなり充分実用的に見える。同じような検証してる人が他にも居て、

Google グループ

この人は精度は1/15ピクセルだと結論している。ま、コンパラブルだわ。

solve-field --verifyの効果

何らかの原因でσが大きい場合、まあ、3ピクセルとかだとパノラマにした時に具合悪いわな、何かできることはないのかと。あるそうです。一旦解いた結果ファイルをもう一回食わせてやるとよいらしい。んで、やってみた。リネームしないと怒られます。

$> cd G
$> mv DSC_4395.wcs DSC_4395_.wcs
$> solve-field --overwrite --verify DSC_4395_.wcs DSC_4395.tif

結果、σは1.608912e+00から8.229623e-01へ改善。へーってなもんです。毎度fv開いてちまちま計算してられないので、Gnu Rのスクリプトをでっち上げた。いや便利だわ。

require (FITSio)

corrfiles <- system ("ls G/DSC_*.corr", intern=T)
sd <- c ()
for (corrfile in corrfiles){
  zz <- file (description = corrfile, open = "rb")
  header0 <- readFITSheader (zz)
  header <- readFITSheader (zz)
  D <- readFITSbintable (zz, header)
  close (zz)
  
  # str (D)
  # str (header)
  # str (parseHdr (header))
  D$colNames
  sd <- c (sd, sd (sqrt *1
}

corrfiles
sd

*1:D$col1 - D$col5)^2 + (D$col2 - D$col6)^2))) summary (sqrt ((D$col\[\[1\]\] - D$col\5)^2 + (D$col\2 - D$col\6)^2