[KEKB Bunch Feedback Group]

PF-AR用DCCT読み出しとビーム寿命/入射レート計算(Japanese)


by とびやま まこと(Makoto Tobiyama)/KEKB ビームモニターグループ

警告
以下の記述に関しては、意図する、しないに関わらず多くの誤り、誤解が含まれていると思われますので、決して信用してはいけません。これを信じて起きた損害に関しては、当方は一切責任を持ちません。


If you need contact with the author, please E-mail makoto.tobiyama@kek.jp.
目次

1.はじめに

今から思えば昔のことですが(とはいっても実は今年の7月までのこと--今は1998年9月15日--勤労意欲がわかないのでこんなものを書いている)、PF-AR(舊號トリスタンAR)で新型ビームDCCTのテストを行っておりました。DCCT回路はAR南制御室にあり、正規のビーム電流値はAR南からAR西へ送り、そこでレベルシフト等の回路及びRF等への分配を行った後、CAMACのA/Dで読むルートとデジボルで読みとるルートに分かれ・・といった複雑な経路をたどり、HIDICやCATV表示に使われていました。この正規のルートの他に、モニター端子を直接ハイブリッドレコーダやデジタルマルチメータで読み、それから電流値、ビーム寿命、入射レートを出すことを試みました。

例によってこの文章を書いている今の時点では、寄る年波には勝てず記憶はきわめて怪しげなので、以下の(これまでも)記述には大嘘が多分に含まれるものと思われます。

DCCT装置からは、ハイインピーダンスのビーム電流信号が出てきています。1V=100mAではなかったかしら(あやしげ)。これをアドバンテスト社デジタルマルチメータR6551(DVM)で読みとり、GP-IBを通してEPICSに渡します。これを直接表示するとともに、シーケンサーでビーム寿命及び入射レートを計算します。

2.Adevantest R6551 DVM用GP-IBデバイスサポート

GP-IBのデバイスサポートですが、難しい機能は一切つかいませんのでカンタンです。以下にGDLで書いたものを示します。
Device(R6551)

%{
}%

NameTable Functions = "DCV", "ACV", "2WR", "4WR", "DCA", "ACA";
NameTable DCVRange = "AUTO", "300mV", "3000mV", "30V", "300V", "1000V";
NameTable ACVRange = "AUTO", "300mV", "3000mV", "30V", "300V", "700V";
NameTable ResisRange = "AUTO", "300 Ohm"," 3000 Ohm", "30 kOhm", "300 kOhm",
                     "3000 kOhm", "30 MOhm", "300 MOhm";
NameTable CurrentRange = "AUTO", "300mA", "3000mA";
NameTable SampleMode  = "Free", "Hold";
NameTable SampleRate  = "Fast", "Middle", "Slow";
NameTable ViewCols    = "3+1/2", "4+1/2", "5+1/2";
NameTable Nulls       = "OFF", "ON";
NameTable Scaling     = "OFF", "ON";
NameTable Filter      = "OFF", "ON";
NameTable AutoZero    = "OFF", "ON", "ONCE";  

EfastTable FunctionsSet = "F1","F2","F3","F4","F5","F6";
EfastTable DCVRangeSet =  "F1,R0","F1,R3","F1,R4","F1,R5","F1,R6","F1,R7";
EfastTable ACVRangeSet =  "F2,R0","F2,R3","F2,R4","F2,R5","F2,R6","F2,R7";
EfastTable TwResisRangeSet =  
     "F3,R0","F3,R3","F3,R4","F3,R5","F3,R6","F3,R7","F3,R8","F3,R9";
EfastTable FwResisRangeSet =  
     "F4,R0","F4,R3","F4,R4","F4,R5","F4,R6","F4,R7","F4,R8","F4,R9";
EfastTable DCCurrentRangeSet = "F5,R0", "F5,R6", "F5,R7";
EfastTable ACCurrentRangeSet = "F6,R0", "F6,R6", "F6,R7";
EfastTable SampleModeSet = "M0", "M1";
EfastTable SampleRateSet = "PR1","PR2","PR3";
EfastTable ViewColsSet = "RE3","RE4","RE5";
EfastTable NullsSet = "NL0","NL1";
EfastTable ScalingSet = "SC0","SC1";
EfastTable FilterSet = "FL0","FL1";
EfastTable AutoZeroSet = "AZ0","AZ1","AZ2";

ParamTable{
        set_function	{rec=mbbo, efast=FunctionsSet, name=Functions}
        set_DCV_range	{rec=mbbo, efast=DCVRangeSet, name=DCVRange}
	set_ACV_range	{rec=mbbo, efast=ACVRangeSet, name=ACVRange}
	set_TwR_range	{rec=mbbo, efast=TwResisRangeSet, name=ResisRange}
	set_Fwr_range	{rec=mbbo, efast=FwResisRangeSet, name=ResisRange}
	set_DCC_range	{rec=mbbo, efast=DCCurrentRangeSet, 
				   name=CurrentRange}
        set_ACC_range   {rec=mbbo, efast=ACCurrentRangeSet,
                                   name=CurrentRange}
        set_Samp_mode   {rec=bo, efast=SampleModeSet, name=SampleMode}
	set_Samp_rate	{rec=mbbo, efast=SampleRateSet, name=SampleRate}
	set_VCols	{rec=mbbo, efast=ViewColsSet, name=ViewCols}
	set_Nulls	{rec=bo, efast=NullsSet, name=Nulls}
	set_Scaling	{rec=bo, efast=ScalingSet, name=Scaling}
	set_Filter	{rec=bo, efast=FilterSet, name=Filter}
	set_Autozero	{rec=mbbo, efast=AutoZeroSet, name=AutoZero}
	get_data	{rec=ai, command="E\r\n",conv="%lf\r\n",leng=511}

}
%{

}%

内容については説明するまでも無いと思います。なお、DCしか使いませんので、上の関数のうち使うのはget_data、set_DCV_range、set_SAMP_rage位です。

3.EPICSデータベース

データベースを作ります。デバイスサポートの関数以外に、ビーム電流への換算とか、ビーム寿命とか、入射レートとかのデータベースが必要ですし、オフセットを除去できるような仕組みを作っておきます。capfastで書いたデータベースを示します。

get_data自身はSCANを「1 second」に設定しましたので、(DVMのサンプルとは非同期に)1秒ごとにデータを読んで来ます。ここで、DVMで読んだ値を100倍すれば、mA単位での電流値になるのですが、オフセットを調整するため、同じ値をaiに入れて置いて、ある時これが「0mA」だと指示したとき、それを下駄にするよう、100*(A-B)とします。ここで、AがDVMで読んだ生の値です。これに加えて、この補正に失敗したとき手で直せるように、calcのCに値をいれることが出来るようにしておきます。結局、calcの計算式は(A-B)*100.0-Cとなります。

このほかの、例えばfit_startとかは、次のビーム寿命等の計算の時に使います。

4.ビーム寿命、入射レート計算用シーケンサ

リングを回っているビームが出す放射光やウエークフィールドが引き起こす脱ガスによるリング真空度の変化が大きくないとすると(この条件はARには全く当てはまらない)、蓄積ビーム電流値は時間とともに以下の式で表されるような形で減少する。

I=I0exp(-t/τ)
この「τ」が(色々な効果を合わせた)ビーム寿命となる。実際にはARの真空系はあまりパワフルとは言えず(控えめな言い方)、ビーム電流の大きいところと小さな所ではかなり真空度が違い、ビーム寿命は上の式からかなり逸脱する(exponentailでなく直線的になる)が、ある区間を限ってみれば、exponentialの式でビーム寿命を定義しても問題はない。

実際には、1秒刻みでサンプルした今から過去に向かっての電流値のデータをいくつか記録しておいて、それをexponentialの式で最小2乗近似ををしてビーム寿命を出すことになります。多くの点数を取ればそれだけ1点毎のノイズ等の影響は小さくなりますが、あまり寿命が短いときに延々と長くフィット点数を取れば、元々の式からの逸脱も大きくなるし、また寿命の急変等にも対応出来なくなります。そこで、「今」のlifeの値をもとに寿命計算のフィット点数を判断する必要があります。また、この点数の切り替えにより(計算が変わるので)計算したlifeに若干の差がでることもあり、下手をすると切り替わりが非常に頻繁に起こってしまいます。そこで、切り替えにはシュミットトリガ的な動作が必要になります。シーケンサーの全ソースコードはここをクリックして下さい。

program ardcct_l2
%{extern double log();}%
%{extern double fabs();
  extern double exp();}%

 float FB;
 float LIFE;
 float RATE;

%{ 
   float ib[500];
   float s_xy,s_x,s_y,s_x2,a1,a2,a3,a4,a5,a6;
   float samp_time = 1.0;
   int   l_10, l_60, l_240;
   int	 llflag1,llflag2,llflag3,llflag4;
}%
logとかabsとかexpはシーケンサーとしてexternalファンクションになりますので、「c言語の」externalだとして宣言してやる必要があります。また、チャンネルと直接お話をしない部分は普通のCの書き方をすればよいので、%{と}%でくくって書きます。
 assign FB to "fb_ar_dcct:dcct";
 assign LIFE to "fb_ar_dcct:life";
 assign RATE to "fb_ar_dcct:rate";
 assign FITSTART to "fb_ar_dcct:fitstart";

 assign P_MON01 to "fb_ar_dcct:dcct";

 monitor P_MON01;

チャンネルへのassignment及びmonitor宣言をします。今回のmonitorはdcct電流値の計算終了です(DVMへのサンプルと同じ)。これは実は値が変わらないと動かないので、本当はちょっとまずい。
ss state_main1
{
	state init
	{
		when()
		{
		printf("init state\n");
		ZERO = 0.0;
		MAXFIT = 500;
		LIFE = 0.001;
		FITSTART = MAXFIT -10;
		FITN = MAXFIT -FITSTART;
		MIN=0.2;
		MIN2=0.003;
		P_MON01P=0;
		samp_time = 1.0;
		llflag1 = 0;
		llflag2 = 0;
		llflag3 = 0;
		llflag4 = 1;
		
		l_10 = 0;
		l_60 = 0;
		l_240 = 0;

		for (i=0; i<200; i++)
		{
			ib[i]=0.01;
		}
		}state sta1
	}

initステートは立ち上がったとき1回だけ動くsteteです。次はsta1ステートへ行き、戻っては来ません。
        state sta1
        {
                when(P_MON01!=P_MON01P)
		{
		s_x = 0.0;
		s_x2 = 0.0;
		s_y = 0.0;
		s_xy = 0.0;
                pvGet(FB);
/*		printf("FB= %f\n",FB); */

		if (FB<MIN)
		 ib[MAXFIT-1]=0.01;
		else 
		 ib[MAXFIT-1] = FB;
最小2乗法で簡単に計算するために、ビーム電流のlogをとります。このため、電流値が0とか負の値は大変困ってしまうので、その時はあきらめて0.01mAに設定してしまいます。
		if ((LIFE > 30.0) && (llflag3))
		   {
		   if ((l_240 > 239) && (!llflag4))
		    {
		     FITSTART=MAXFIT-240;
		     pvPut(FITSTART);
		     llflag4=1;
		     l_240 = 240;
		     }
		   else
		    {
		    l_240++;
		    }
		   }
		else if (LIFE > 15.0)
		   {
		   if ((l_60 > 59) && (!llflag3) || (llflag4))
		      {
		      FITSTART=MAXFIT-60;
		      pvPut(FITSTART);
		      llflag3 = 1;
		      llflag4 = 0;
		      l_60 = 60;
		      l_240 = 60;
		      }
		   else 
		     {
		     l_60++;
		     }
		   }
		else if (LIFE> 1.0)
		   {
		   llflag3 = 0;
		   llflag4 = 0;
		   l_60 = 0;
		   l_240 = 0;

		   FITSTART=MAXFIT-10;
		   pvPut(FITSTART);
		   }
		else if (LIFE> 0.0)
		   {
		   llflag3 = 0;
		   llflag4 = 0;
		   FITSTART = MAXFIT - 5;
		   pvPut(FITSTART);
		   l_60 = 0;
		   l_240 = 0;
		   }
		else
		  {
		  llflag3 = 0;
		  llflag4 = 0;
		  l_60 = 0;
		  l_240 = 0;
		  FITSTART = MAXFIT - 3;
		  pvPut(FITSTART);
		  }  		   

今のビーム寿命をもとにフィットする点数を決めます。基本的に、寿命が30分以上なら240点(4分)、15分以上なら60点(1分)、1分以上なら10点、それ以下でビームがふえていなければ5点、ビームがふえているなら(入射中のはず)3点としました。また、一気にモード変更が起きると拙いので、その寿命が何回か続くと次のモードへ移行するようにしました。ただし、lifeが急に長くなる理由はあまりないと思いますが、不幸にしてlifeが急に短くなる理由はいくらでもありますので、短くなったときは問答無用で短寿命モードへ移行するようにしました。今回のARの運転では、EPICSデータベース上では今がどんな運転モードか(入射中か、2.5GeV蓄積中か、加速中か、6.5GeV蓄積中か)分からないので、やはり入射近辺ではライフ計算モード移行がうまくいかず、しばらくがたがたが続くこともありました。この辺はなかなか奥が深そうで、下手をするとビーム寿命計算がライフワークになってしまうかも知れません(そいつはご勘弁を)。
		for (i=FITSTART;i<MAXFIT; i++)
		{
			s_x  += i*samp_time;
			s_x2 += (i*samp_time)*(i*samp_time);
			s_y  += log(ib[i]);
			s_xy += (i*samp_time)*log(ib[i]);
		}
		FITN=MAXFIT-FITSTART;
		a1 = FITN*s_x2 - s_x*s_x;
		a2 = FITN*s_xy - s_x*s_y;
		a3 = a2/a1;

		RATE = a3*ib[FITSTART];
		pvPut(RATE);

		a3 = -a3*60;

/*		printf(" a3 and ABS =%f  %f\n",a3,fabs(a3));  */
		

		if (fabs(a3)<MIN2)
		   a4 = 0;
		else 		   
		   a4 =1.0/a3;
		LIFE = a4;
		pvPut(LIFE);
		P_MON01P=P_MON01;
		for (i=0;i<MAXFIT-1; i++)
		{
			ib[i] = ib[i+1];
		}
                }state sta1
        }
}
実際の最小2乗計算をしているところです。線形フィットをして、傾き分の1が寿命になります。これを分の単位に換算して終わりです。

5.MEDMによる表示

私は未熟者でmedmしか使えませんので、medmでbeam current、life、injection rateを表示するパネルを作りました。しかしながら、medmはお手軽ではありますが融通が利きませんので、一般に不評でした。皆様から寄せられましたご批判は私も同様に感じている所でして、次からは多分ぱいそんを使ったパネルになるものと(少しだけ)ご期待下さい。
動作中の画面を記録しておくのを忘れていましたので、設計画面のみ示します。
medm image

6.おわりに

PF-ARの3月立ち上げ用に作ったDCCT読みとりテストセットについて説明しました。シーケンサをまともに動作させるまでに、コントロールグループの山本氏、中村氏に色々ご助言を賜りました。感謝いたします。

実際には、運転中にシーケンサを止め、改造し、再ロードしたくなることが多くあると思われますが、不思議なことに今まで再ロードに成功した試しがありません。よくいって動作しないか、あるいはIOCがハングアップしてしまうといった状態です。まあ、あきらめてリブートする方が早いかも知れませんが、運転中にDCCTの情報が無くなるといった事態が起きて良い物か、どうでしょう。


Makoto Tobiyama
15/Sep/98

Return to FB Home Page...