|
SPECの生データをIgorProで表示，(x,y,err)形式のデータとしても保存
(1) rawdataファイルの生成 SPECの生データをスキャン番号別のファイル***.txtとして書き出す．
以下の内容のファイルspec2txtを作り，**.txtを作りたいrawdataディレクトリに置く．
- BEGIN{FLAG=0}
$1=="#S" {FLAG=1; NAME=$2; printf("") > (NAME".txt")};
$0=="" {FLAG=0};
{if(FLAG==1) print $0 >> (NAME".txt")}
{if(FLAG==0) close(NAME".txt")}
rawdataディレクトリで次を実行．生データは1つ上にあるとする．
- >awk -f spec2txt ../150928
(2) rawdataファイルをIgorProで読み込んでデータを表示する
- Graph0での操作と Table1, Table2
- Graph0のL_Pathに***.txtファイルが入ったrawdataディレクトリのパスを入力．
- Scan No. の数字を変えると対応するスキャン番号の結果がGraph0に表示される．
- ***.txtファイルから読み取ったスキャン情報はTable2に表示される．
- データはTable 1に表示されるので，それを見て，X軸，Y軸にしたいカラム番号をそれぞれX-axis（1列目dat[][0]のカラム番号は0である），Y-axisに指定する．
 
(3) rawdataファイルをIgorProで読み込んでデータを(x,y,err)形式にして保存する IgorProは上と同じViewScansを使用．C-PLOTでやっていることをIgorProでもできるようにした．
- Graph1とTable0での操作
- Graph1のS_Pathに(x,y,err)形式にした***.datファイルを保存するためのdataディレクトリのパスを入力
- Table0のscn[][0]列にスキャン番号を記入．複数のスキャンを重ね合わせて1つにするときは，scn[][1]~scn[][5]列にも順に記入．ただし，最後は0（3つのスキャンを重ね合わせるときは，scn[][0]~scn[][2]に3つのスキャン番号を記入し，scn[][3]は0にしておく）．scn[][5]列までで足りないときはscn[][6]列以降を付け足す（redimension/N=(行数,列数) scn）．
- xの差がbindxに指定した値以下であれば同じxとみなして重ね合わせる．
- x軸（横軸の変数）にする列番号をxcol，y軸（データ）にする列番号をdatcol，モニター値の列番号をmoncol，1秒間のモニター値（統一する値）をmon1sに記入．moncolに-2と書くと，後ろから2列目の意味．datcolも同様．Table1のdata[][0]をx軸にする場合，1列目なので，xcolは1にする（0にしない）．
- (x,y,err)形式の***.datファイルを作るときは，Graph0のi1とi2に，それぞれTable0で実行したい最初の行番号と最後の行番号を入れ，calcを押す．すると，S_Pathに指定したdatディレクトリにScan***.dat（***にスキャン番号）の名前で(x,y,err)形式になったデータファイルが保存される．
 
データ処理のためのIgorPro関数
C-PLOTまたはIgorProで作られた大量のscan***.datファイルから必要なものを選んで読み込んできて，プロットしたりフィットしたり，いろいろな解析をするための関数をまとめておく．以下の関数をProcedureに埋め込み，コマンドラインで使う．データが大量なので，以下のような作業を一つ一つメニューからコマンドを選んでやると膨大な時間がかかる．いかに関数やマクロを組んでまとめて処理するかを考えて作業することが肝要． 注：以下の Function . . End の部分をIgorProのProcedureにコピーペーストすればいいのだが，ここはWeb上のテキストであるため，Tabや連続スペースがそのままうまく反映されない．普通は，区切りごとに段差を設けて見やすい形にするのだが，ここは全部1カラム目から書いてあって見にくい．IgorProにペーストしたあと，Tabやスペースを入れて見やすくしたほうがいい． (1) スキャン番号waveの作成 一連のスキャン番号が入力されたwaveをData -> Make Wavesで作っておく．例では，名前がscantth（#302–350)とscanth(#353–365)としてある．測定のまとまりごとにうまく名前を使いわけるといい．この例ではtthスキャンだけとthスキャンだけをまとめたものになっている．POLという変数がパラメータになっていて，POL=-90から+90まで15きざみでのtthスキャンとthスキャンの結果である．点数は13点であるが，点数を数えるのは面倒なので，適当に20点くらいで作っておき，余分なところは後で消すほうが楽．
•Make/N=20/D scantth,scanth
•edit scantth,scanth
scantthは#302から4ずつ，scanthは#353から1ずつ増えているので，次のようにすると一気にできる．
•scantth=302+4*p; scanth=353+p
(2) データを一気に読み込む Function load(snums,path)
wave snums
String path
String wname=NameOfWave(snums)+"_"
variable i=0
do
LoadWave/G/D/A/Q path+"scan"+num2str(snums(i))+".dat"
// LoadWave/G/D/A path+"scan"+num2str(snums(i))
duplicate/O wave0 xx
duplicate/O wave1 yy
duplicate/O wave2 err
duplicate/O xx $wname+num2str(i)+"X"
duplicate/O yy $wname+num2str(i)+"Y"
duplicate/O err $wname+num2str(i)+"E"
killwaves wave0,wave1,wave2
i+=1
while (i < numpnts(snums))
End
使い方：
•load(scantth,"Macintosh HD:Users:tmatsu:Desktop:example:data:")
scantthはスキャン番号が入ったwave名，”Mac---data:”の部分はscan***.datの名でデータ（C-PLOTの出力）が保存されているファイルパス．ファイルパスの書き方はDataメニューのLoad Wavesで何かファイルを選択してみるとわかる．scan***.datはエラーバー付きの(x, y, s)を一組とする3列データ．
新しくデータが更新されたときなど，データを1個だけ読み込むときは，次の関数を使うといい．
Function load1(snums,i,path)
wave snums
variable i
String path
String wname=NameOfWave(snums)+"_"
LoadWave/G/D/A path+"scan"+num2str(snums(i))+".dat"
// LoadWave/G/D/A path+"scan"+num2str(snums(i))
duplicate/O wave0 xx
duplicate/O wave1 yy
duplicate/O wave2 err
duplicate/O xx $wname+num2str(i)+"X"
duplicate/O yy $wname+num2str(i)+"Y"
duplicate/O err $wname+num2str(i)+"E"
killwaves wave0,wave1,wave2
End
使い方：
•load1(scantth,12,"Macintosh HD:Users:tmatsu:Desktop:example:data:")
scantthのp=12のデータ（scan350.dat）だけが読み込まれる．
(3) すべてのスキャンを1枚のグラフにプロットする Function dispall(snums)
wave snums
String wname=NameOfWave(snums)+"_"
variable i=1
Display $wname+num2str(0)+"Y" vs $wname+num2str(0)+"X"
ErrorBars/T=0.5/L=0.5/Y=2 $wname+num2str(0)+"Y" Y,wave=($wname+num2str(0)+"E",$wname+num2str(0)+"E")
do
Appendtograph $wname+num2str(i)+"Y" vs $wname+num2str(i)+"X"
ErrorBars/T=0.5/L=0.5/Y=2 $wname+num2str(i)+"Y" Y,wave=($wname+num2str(i)+"E",$wname+num2str(i)+"E")
i+=1
while (i < numpnts(snums))
ModifyGraph tick=2,mirror=1,standoff=0,width=283,height=170
End
使い方：
•dispall(scantth)
scantthに含まれる#302から#350までのすべてのデータが1つのウィンドウにプロットされる．
(4) 1つのスキャンだけプロットする． Function graph(snums,i)
wave snums
variable i
String wname=NameOfWave(snums)+"_"
Display $wname+num2str(i)+"Y" vs $wname+num2str(i)+"X"
ErrorBars/T=0.5/L=0.5/Y=2 $wname+num2str(i)+"Y" Y,wave=($wname+num2str(i)+"E",$wname+num2str(i)+"E")
ModifyGraph tick=2,mirror=1,minor=1,sep=10,standoff=0
ModifyGraph mode=4,marker=19,msize=3,opaque=1
ModifyGraph lsize=0.5,rgb=(1,3,39321)
End
使い方：
•graph(scantth,5)
scantthのp=5のデータ（scan322.dat）をプロットする．エラーバー付き．
(5) 別のスキャンを追加でプロット（Append）する Function app(snums,i)
wave snums
variable i
String wname=NameOfWave(snums)+"_"
Appendtograph $wname+num2str(i)+"Y" vs $wname+num2str(i)+"X"
ErrorBars/T=0.5/L=0.5/Y=2 $wname+num2str(i)+"Y" Y,wave=($wname+num2str(i)+"E",$wname+num2str(i)+"E")
ModifyGraph mode($wname+num2str(i)+"Y")=4,msize($wname+num2str(i)+"Y")=3,opaque($wname+num2str(i)+"Y")=1,lsize($wname+num2str(i)+"Y")=0.5
End
使い方：Appendしたいウィンドウを最前面に置いてから次を実行する．
•app(scantth,8)
scantthのp=8のデータ（scan334.dat）をAppendする．エラーバー付き．
(6) すべてのスキャンを別々のウィンドウにプロットし，別々のコメントをつける すべてのスキャンを別々のウィンドウでプロットし，ウィンドウごとに測定パラメータを変えたコメントを表示する．エラーバー付き．測定結果を全部プリントアウトして，一目で全スキャンの状況が見られるようにするためのもの．このプリントアウトは，データの解析，実験レポートの作成，プレゼン準備，論文作成など，後々いろいろな場面で参照することが多い．測定の異常に気がつくためにも必要．記入ミスのないよう，ちゃんと作っておく．
Function allgraph(snums,param)
wave snums, param
String wname=NameOfWave(snums)+"_"
variable i=0
do
graph(snums,i)
Label left "Intensity (cnts/s)";Label bottom "two-theta (deg.)"
SetAxis left 0, 100000 ;SetAxis bottom 84,89;
ModifyGraph margin(left)=51,margin(bottom)=43,margin(top)=14,margin(right)=14
ModifyGraph width=184,height=136,lsize=0.5,rgb=(0,0,0), gfSize=12
String text1="(3 0 1.5)\rT=10 K\r"
String text2="\\F'Symbol's-p\\F'Helvetica''"
String text3="\rPol="+num2str(param(i))+" °"
TextBox/C/N=text0/F=0/A=LT/X=5/Y=5/B=1 text1+text2+text3
i+=1
while (i < numpnts(snums))
End
使い方：p=0~12までのスキャンがpol=-90から90まで15きざみのスキャンに対応しているので，パラメータであるpolの値が入ったwaveを作る．duplicateするのがいい．
•duplicate scantth pol
•pol=-90+p*15; edit pol
上の関数のLabel left . . .からTextBox . . . の行までが各グラフのスタイルやコメント内容を指定する部分になっているので，ここを好きなように書き換える．String text3= . . . が，グラフごとにコメントが変わる部分．試料名，エネルギー，磁場，スリット，アッテネータなどの情報も記入したいが，あまりにごちゃごちゃしてくるので，ここは最低限重要なものにとどめ，共通パラメータはレイアウトに記入するか，プリントアウトに手書きで記入してもよいだろう．
•allgraph(scantth,pol)
1回で満足のいくプロットができて成功することはあまりない．納得がいくまで何回か修正する．大量にできてしまったグラフウィンドウを消去するには，Windows -> Control -> Window BrowserでWindow Browserを開き，Graphsで消去したいグラフを選び，Killする．あるいはCommand+Option+Wで消去する．
1枚の紙にまとめてプリントアウトするには，メニューのNew LayoutでTileにチェックを入れ，まとめて表示したいグラフを選んで実行するとよい．A4用紙1枚に5行3列で15枚のグラフまでなら，なんとか読むに支障がない程度に表示できる．LayoutのPage Setupで縮小度を70%くらいにするといい．
(7) 1つのスキャンについてy値とs値（エラーバー）を定数倍する Function mult(snums,fac,i)
wave snums
variable fac,i
String wname=NameOfWave(snums)+"_"
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
yy=yy*fac
err=err*fac
duplicate/O xx $wname+num2str(i)+"X"
duplicate/O yy $wname+num2str(i)+"Y"
duplicate/O err $wname+num2str(i)+"E"
End
使い方：
•mult(scantth,0.1,0)
scantthのp=0のデータ（scan302.dat）のy値とs値が0.1倍される．
(8) すべてのスキャンについてy値とs値（エラーバー）を定数倍する Function multall(snums,fac)
wave snums
variable fac
variable i=0
do
mult(snums,fac,i)
i+=1
while(i<numpnts(snums))
End
使い方：
•multall(scantth,0.1)
scantthの#302から#350までのすべてのデータのy値とs値が0.1倍される．
(9) 1つのスキャンについて，xの関数をyとsにかける たとえば，エネルギースキャンの吸収補正がこれにあたる．
Function multmu(snums,i,muE,mu)
wave snums
wave muE,mu
variable i
String wname=NameOfWave(snums)+"_"
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
yy=yy*interp(xx,muE,mu);
err=err*interp(xx,muE,mu);
duplicate/O xx $wname+num2str(i)+"X"
duplicate/O yy $wname+num2str(i)+"Y"
duplicate/O err $wname+num2str(i)+"E"
End
使い方：escanという名のwaveのp=3のデータがエネルギースキャンになっているとし，それについて，xwave=muE, ywave=muで表される吸収係数をかけるとき，次のようにする．muEとmuはあらかじめ別に作っておく必要がある．
•multmu(escan,3,muE,mu)
(10) 1つのスキャンを yy vs xx Graph にプロットする まとめてフィッティングをする前の試しフィットをするための作業．
Function dispxy(snums,i)
wave snums
variable i
String wname=NameOfWave(snums)+"_"
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
End
使い方：まず，yy vs xx Graphを表示しておく．xx, yy, errはloadでも使われており，いろいろなところで汎用的に使うwaveである．
•display yy vs xx; ErrorBars/T=0.5/L=0.5/Y=2 yy Y,wave=(err,err)
次に，表示したいスキャンを指定してdispxyを実行する．
•dispxy(scanth,0)
yy vs xx Graphにscanthのp=0のデータ（scan353.dat）が表示されるので，これについて，適当なフィッティング関数を決め，試しにフィットしてみる．うまくいけば，その関数を使って，次の(11)ですべてのデータをまとめてフィットする．
(11) すべてのスキャンについて一括フィッティングを行い，パラメータの変化をプロットする Function allfit(snums,W_coef,W_sigma)
wave snums, W_coef,W_sigma
variable nparam
variable i=0
String wname=NameOfWave(snums)+"_"
nparam=numpnts(W_coef)
Make/N=(nparam,numpnts(snums))/D/O W_par
Make/N=(nparam,numpnts(snums))/D/O W_sig
do
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
//print i
//FuncFit/H="0000000" Gausfit W_coef yy /X=xx /W=err /I=1 /D
CurveFit/H="1000"/NTHR=0/TBOX=0 lor yy /X=xx /W=err /I=1 /D
W_par[][i]=W_coef[p]
W_sig[][i]=W_sigma[p]
i+=1
while (i < numpnts(snums))
duplicate/O W_par $wname+"Wpar"
duplicate/O W_sig $wname+"Wsig"
End
使い方：先に(10)で試しフィットをしてフィッティングの方針を決め，次にここで連続フィットをする．W_coefはフィッティングパラメータが入るwave名，W_sigmaは標準偏差sigmaが入るwave名．CurveFIt . . . はIgorPro組み込み関数を用いたフィット．FuncFit . . . はユーザー定義関数によるフィット．この部分を自分がやりたいように書き換える．上の例では，IgorPro組み込みのLorentzianでW_coef[0]（バックグラウンド）を固定してフィットする指定になっている．1つの行で//から後はコメント扱い．いろいろなフィッティングパターンを記入しておくといい．次を実行する．
•allfit(scanth, W_coef, W_sigma)
W_coefとW_sigmaは(10)の試しフィットで既に作られており，それをここでも使用する．これを実行すると，scanthのデータが次々とyy vs xx Graphに上書きされ，そこでフィッティングが行われていく．yy vs xx Graphを最前面に出してから実行し，フィッティングがうまくいっているかどうか，様子を観察する．結果のパラメータはscanth_Wparに，標準偏差はscanth_Wsigに書き込まれている．
(12) フィッティングパラメータの変化をプロットする 中心値の変化をプロットしてみよう．まず，フィットで得られた中心値を入れるためのwaveを作り，パラメータであるpolと同じ表に表示する．polの表を最前面に出して次を実行する．•duplicate pol centh; appendtotable centh
中心値はp=2の欄に入っているので，次のように実行する．
•centh=scanth_Wpar[2][p]
これをプロットする．
•display centh vs pol
(13) フィッティング結果を用いて積分強度を求め，強度変化をプロットする Function background(snums,nbg)
wave snums
variable nbg //number of points at each ends of the data taken for the background estimation
variable i=0
String wname=NameOfWave(snums)+"_"
Make/N=(numpnts(snums))/D/O bg,bger
do
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
duplicate/O err err2; err2=err^2
//nbg=3*(i==0)+5*(i>0)*(i<=6)+8*(i>6)
bg[i]=(sum(yy,0,nbg-1)+sum(yy,numpnts(yy)-nbg,numpnts(yy)-1))/2/nbg
bger[i]=sqrt((sum(err2,0,nbg-1)+sum(err2,numpnts(yy)-nbg,numpnts(yy)-1)))/2/nbg
i+=1
while (i < numpnts(snums))
//duplicate/O $wname+"Wpar" W_par; duplicate/O $wname+"Wsig" W_sig
//bg=W_par[0][p]; bger=W_sig[0][p]
duplicate/O bg $wname+"BG"
duplicate/O bger $wname+"BGerr"
End
Function integsum(snums)
wave snums
String wname=NameOfWave(snums)+"_"
variable i=0
Make/N=(numpnts(snums))/D/O interr,int
do
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
duplicate/O $wname+"BG" bg
duplicate/O $wname+"BGerr" bger
variable intsum=0, intsum_err=0, ii=0
do
intsum=intsum+(yy[ii]+yy[ii+1])*(xx[ii+1]-xx[ii])/2
intsum_err=sqrt(intsum_err^2+(err[ii]^2+err[ii+1]^2)*(xx[ii+1]-xx[ii])^2/2)
ii+=1
while(ii<numpnts(xx)-1)
int[i]=intsum-bg[i]*(xx[numpnts(xx)-1]-xx[0])
interr[i]=sqrt(intsum_err^2+(bger[i]*(xx[numpnts(xx)-1]-xx[0]))^2)
i+=1
while (i < numpnts(snums))
duplicate/O int $wname+"INT"
duplicate/O interr $wname+"INTerr"
End
Function integfit(snums)
wave snums
String wname=NameOfWave(snums)+"_"
variable i=0
Make/N=(numpnts(snums))/D/O interr,int,FWHMerr,FWHM
duplicate/O $wname+"Wpar" W_par
duplicate/O $wname+"Wsig" W_sig
// for Igor built in Gaussian
// int =W_par[1][p]*abs(W_par[3][p])*sqrt(pi)
// interr=sqrt(W_par[3][p]^2*W_sig[1][p]^2+W_par[1][p]^2*W_sig[3][p]^2)*sqrt(pi)
// FWHM=2*abs(W_par[3][p])*sqrt(ln(2))
// FWHMerr=2*abs(W_sig[3][p])*sqrt(ln(2))
// for Igor built in Lorentzian
int =W_par[1][p]/sqrt(abs(W_par[3][p]))*pi
interr=sqrt(W_sig[1][p]^2+W_par[1][p]^2*W_sig[3][p]^2/W_par[3][p]^2/4)/sqrt(abs(W_par[3][p]))*pi
FWHM=2*sqrt(abs(W_par[3][p]))
FWHMerr=abs(W_sig[3][p])/2/sqrt(abs(W_par[3][p]))
duplicate/O int $wname+"INT"
duplicate/O interr $wname+"INTerr"
duplicate/O FWHM $wname+"FWHM"
duplicate/O FWHMerr $wname+"FWHMerr"
End
使い方：全データの台形積分を積分強度として扱う場合はintegsumを使う．まず， 関数backgroundを実行し，データの両サイド何点かずつからバックグラウンドを見積る（結果はbgとbgerに入る）．
•background(scanth,3)
•integsum(scanth)
実行すると，scanth_INTとscanth_INTerrという名前のwaveができているので，POLと同じ表に列を追加して表示する．また，横軸をpol，縦軸を積分強度にとってプロットする．
•appendtotable scanth_INT, scanth_INTerr
•display scanth_INT vs pol; ErrorBars/T=0.5/L=0.5/Y=2 scanth_INT Y,wave=(scanth_INTerr,scanth_INTerr)
Lorentzianでのフィットがきれいにできているときは，フィッティング結果のパラメータから積分強度が計算できる．IgorPro組み込みのLorentizanの場合の積分強度とその誤差の評価式が上の関数integfitの中に書かれているので，それを使って計算する．結果は台形積分と同様にscanth_INTとscanth_INTerrに上書きされる．
•integfit(scanth)
ユーザー定義関数によるフィットで，パラメータがそのまま積分強度になっているときは，たとえばパラメータwaveのp=1のところが積分強度だとすると，次のようにすればよい．
•scanth_INT=scanth_Wpar[1][p]; scanth_INTerr=scanth_Wsig[1][p];
(14) 積分強度の規格化 あるスキャンシリーズの積分強度を別のスキャンシリーズの積分強度で割って規格化する．たとえば(100)反射のtheta-scanだけ集めたscanth1というwaveと，(200)反射のtheta-scanだけ集めたscanth2というwaveがあって，(100)の積分強度を(200)の積分強度で割りたい場合などがこれにあたる．
Function normint(snums1,snums2)
wave snums1,snums2
String wname1=NameOfWave(snums1)+"_"
String wname2=NameOfWave(snums2)+"_"
duplicate/O $wname1+"INT" int
duplicate/O $wname1+"INTerr" interr
duplicate/O $wname2+"INT" intF
duplicate/O $wname2+"INTerr" intFerr
duplicate/O $wname1+"INT" NInt
duplicate/O $wname1+"INTerr" NInterr
NInt=int/intF
NInterr=sqrt(interr^2+int^2/intF^2*intFerr^2)/intF
duplicate/O NInt $wname1+"NInt"
duplicate/O NInterr $wname1+"NIerr"
End
使い方：(12)で2つのスキャンシリーズすべての積分強度とその誤差の計算が終わってから次を実行する．もちろん，scanth1とscanth2の点数は一致していなければならない．
•normint(scanth1,scanth2)
実行後，規格化された(100)の積分強度が入ったscanth1_NIntと，その誤差が入ったscanth1_NInterrというwaveができている．
(15) ピークトップを1に規格化 スキャンの最大値がすべて1になるよう規格化する．最大値は_maxがついたwaveに出力される．
Function normto1(snums)
wave snums
Variable i=0
String wname=NameOfWave(snums)+"_"
Make/D/N=(numpnts(snums))/O W_maxs
do
duplicate/O $wname+num2str(i)+"X" xx
duplicate/O $wname+num2str(i)+"Y" yy
duplicate/O $wname+num2str(i)+"E" err
W_maxs[i]=wavemax(yy)
yy=yy/W_maxs[i]
err=err/W_maxs[i]
duplicate/O yy $wname+num2str(i)+"Y"
duplicate/O err $wname+num2str(i)+"E"
duplicate/O W_maxs $wname+"_max"
i+=1
while(i<numpnts(snums))
End
(16) bindata : データを横軸方向に何点かまとめ，平均化する Function bindata(xwave,ywave,errwave,nbin)
wave xwave,ywave,errwave
variable nbin
duplicate/O xwave xx
duplicate/O ywave yy
duplicate/O errwave err
duplicate/O err err2; err2=err^2
Make/O/N=(numpnts(xx)/nbin)/D xbin,ybin,errbin
variable ii=0
do
xbin[ii]=sum(xx,ii*nbin,(ii+1)*nbin-1)/nbin
ybin[ii]=sum(yy,ii*nbin,(ii+1)*nbin-1)/nbin
errbin[ii]=sqrt(sum(err2,ii*nbin,(ii+1)*nbin-1))/nbin
ii+=1
while(ii*nbin<=numpnts(xx)-nbin)
duplicate/O xbin xx
duplicate/O ybin yy
duplicate/O errbin err
End
(17) azifit : 基本反射のアジマス角依存性フィット Function azifit(w,x)
Wave w; Variable x
Variable r
if(x<250)
r=w[0]+w[1]*sin(2*pi*x/360+w[2])
else
r=w[3]+w[1]*sin(2*pi*(x-360)/360+w[2]-pi/2)
endif
return r
End |
