LB2045の校正

LB2045を日本分析化学会の放射能分析用玄米認証標準物質(写真:外観容器に詰めたところ測定のようす認証書)で校正した。

出力

認証書によれば 137Cs は (210±13)Bq/kg,134Cs は (141±9)Bq/kg,40K は (75±7)Bq/kg で,セシウムは小さめに,カリウムは大きめに出ている。

なお,認証書によれば2012年6月1日の値であるから,今日(10月5日)は126日を経ており,正確には次のようになり,134Cs : 137Cs 比は0.603ほどである。

210 * 0.5^(126/(365.2422*30.08))
# [1] 208.3372
141 * 0.5^(126/(365.2422*2.0652))
# [1] 125.584
125.584 / 208.3372
# [1] 0.602792

まずRでデータを読み込む:

sp = scan("http://okumuralab.org/~okumura/stat/data/121005rice.txt",
          skip=41, nlines=1024, comment.char=";")

測定時間は10755秒なので,10755を掛けて四捨五入してカウント値に直す:

sp = floor(10755 * sp + 0.5)

これをオレンジでプロットする:

ch = 0:1023          # チャンネル番号
par(mgp=c(2,0.8,0))  # Rのグラフィックオプションを好みでいじる
plot(ch, sp, type="l", xlab="Channel", ylab="Counts", xlim=c(0,400), col="#f39800")

ついでにバックグラウンドも読み込み,青で重ね書きする:

bg = scan("http://okumuralab.org/~okumura/stat/data/121005BaseData.txt",
          skip=36, nlines=1024, comment.char=";")
bg = floor(bg * 10801 + 0.5)
points(ch, bg*(10755/10801), type="l", col="#0068b7")
玄米認証標準物質のスペクトル

カリウムはほとんどバックグラウンドのものであることがわかる。

これをRでフィットする:

library(nlme)
data = data.frame(ch=ch, sp=sp)
r = 260:450
fit = gnls(sp ~ (c + d*(ch-350) + b*(ch-350)^2
                   + (e1*(0.08338/0.9226)*dnorm((ch-0.9314*m1)/(0.9226*s)) +
                      e1*(0.15373/0.9276)*dnorm((ch-0.9415*m1)/(0.9276*s)) +
                      e1*(0.9762/0.956)*dnorm((ch-m1)/(0.956*s)) +
                      e2*0.851*dnorm((ch-m2)/s) +
                      e1*(0.8546/1.0967)*dnorm((ch-m3)/(1.0967*s)) +
                      e1*(0.08688/1.101)*dnorm((ch-1.00765*m3)/(1.101*s))) / s)
                * (ch/350)^u,
           data=data, subset=r,
           start=list(c=200,d=-3,b=0.02, # 初期値
                      e1=10000,e2=10000, # 収束しないときは適宜直す
                      m1=313,m2=342,m3=407,s=9,u=-1.2),
           weights=varPower(fixed=0.5),
           control=list(nlsTol=1e-4)) # 収束しないときは1e-3に直す
a = coef(fit)
plot(ch, sp, type="l", xlab="Channel", ylab="Counts", xlim=c(0,500), col="#f39800")
points(ch[r], fitted(fit), type="l")
points(ch[r], type="l",
       (a['c']+a['d']*(ch[r]-350)+a['b']*(ch[r]-350)^2)
       * (ch[r]/350)^a['u'])
a
#             c             d             b            e1            e2 
#  9.583886e+01 -3.852270e-01  3.052728e-03  9.399035e+03  1.631090e+04 
#            m1            m2            m3             s             u 
#  3.161353e+02  3.442363e+02  4.100476e+02  8.995466e+00 -1.357540e+00 
玄米認証標準物質のスペクトルのフィット

肝心のcalibration factorは,

(210 * 0.377) / (sum(a['e2']*0.851*dnorm((ch[r]-a['m2'])/a['s'])/a['s']*(ch[r]/350)^a['u']) / 10755)
# [1] 59.90988

つまり 59.9 Bq/cps である。LB2045 のページで書いたように,35 Bq/cps あるいは 47 Bq/cps という値を想定していたが,それよりかなり大きめの値になった。

なお,134Cs : 137Cs 比は

a['e1']/a['e2']
#        e1 
# 0.5762427 

となり,最初に計算した値 0.603 より若干小さくなった。


Last modified: 2012-10-05 20:41:08