Lumericalの複素データを2π分の実部データに変換

Lumericalのfrequency monitorで取得して複素数二次元データの形で出力した電場分布をIgorデフォルトの実部表示で見るとある一瞬の位相しか見れないので指定した分割数で2π分の位相を確認するためのProcedureです。
やっていることは大体Lumericalの公式ページ(https://optics.ansys.com/hc/en-us/articles/360034915793-Making-a-CW-movie-from-a-frequency-monitor)がMatlabでやっている手法のmovie出力などを省いたものなのでMatlabが使える人はそちらを使ったほうが良いです。

出力した全ての波長に対して全位相を出力しようとするとデータ数が膨大になってしまうため、あらかじめ何セル目の周波数のデータを見たいか決めて指定する手法をとっています。
今回は全電場分布及び強度平均も見たかったので長々と書いてしまっていますが基本的には

for(i = 0; i < Num_P; i += 1)
pEx[][][i] = real(exp(sqrt(-1)*(2*pi/Num_P*i))*NormalEx[p][q][Num_wl])
endfor

の場所だけしか重要な点は無いので一方向の電場の特定の位相が見たい場合は出力された複素データに対してe^(i*(phase))を乗算して実部をとってしまえば大丈夫だと思います。
X軸やY軸とセルの数が同じだとIgorは画像表示が面倒なため最後にRedimensionでセル数を1減らしてしまっています。セル数が多い場合にはおおむね問題ありませんが厳密には正しくない処理なので論文等に使用する倍は軸の値をちゃんと書き換える必要があります。

Function MakePhaseMatrix()

Variable nm
String	obj = "antenna_L"+num2str(nm)+"nm"	//Lumericalで出力時につける共通の名前

Variable i
Variable Num_P =16			//2πの位相を何分割するか
Variable Num_wl =29			//見たい周波数がFaxisの何cell目か

String Folderpath,FilePath,Fullpath,loadName	//フォルダ名、ファイル名指定、周波数軸読み込み、波長軸変換
	Folderpath	= "D:PD1:Lumerical:Exported:chiral_SPP:"
	FilePath		= "Faxis"
	Fullpath		= FolderPath+FilePath+".mat"
	loadName = "Faxis"
		MLLoadWave /O/M=2/Y=4/S=2/N=$(LoadName) FullPath
		Wave FWave = $(loadName)
		Variable Num_f = rightX(FWave	)
		Duplicate/o FWave WLaxis 	
			WLaxis = 299792458/FWave

String XName,YName,Xload,Yload	//X軸、Y軸読み込み
	XName = "Xaxis"+obj
	YName = "Yaxis"+obj
	Xload = "Xaxis"
	Yload = "Yaxis"
		Fullpath		= FolderPath+XName+".mat"
			MLLoadWave /O/M=2/Y=4/S=2/N=$(Xload) FullPath
		Fullpath		= FolderPath+YName+".mat"
			MLLoadWave /O/M=2/Y=4/S=2/N=$(Yload) FullPath
		Wave Xaxis = $(Xload)
		Wave Yaxis = $(Yload)
Variable Num_x,Num_y 
	Num_X = rightX(Xaxis)		
	Num_Y = rightX(Yaxis)		

String magplane				
	magplane = "mag_"+obj
		Make/o/d/n = (Num_X,Num_Y,Num_F) $(magplane)
		Wave magnitude = $(magplane)
			magnitude = 0
		Duplicate/o magnitude re_f,im_f
		magnitude = 0

String phase_Ex,phase_Ey,phase_Ez,magphase
	phase_Ex = "ExPhase_"+obj
	phase_Ey = "EyPhase_"+obj
	phase_Ez = "EzPhase_"+obj
	magphase = "magPhase_"+obj
		Make/o/d/n = (Num_X,Num_Y,Num_P) $(phase_Ex), $(phase_Ey), $(phase_Ez), $(magPhase)
		Wave pEx = $(phase_Ex)
		Wave pEy = $(phase_Ey)
		Wave pEz = $(phase_Ez)
		Wave pMag = $(magphase)
				
String N_Ex,N_Ey,N_Ez,N_Hx,N_Hy,N_Hz
	N_Ex = "Ex_"+obj
	N_Ey = "Ey_"+obj
	N_Ez = "Ez_"+obj

	Fullpath		= FolderPath+N_Ex+".mat"
	loadName = "N_"+N_Ex
		MLLoadWave /O/M=2/Y=4/S=2/N=$(LoadName) FullPath
		Wave/C NormalEx = $(loadName)
			re_f = real(NormalEx)
			im_f = imag(NormalEx)
			magnitude += re_f^2+im_f^2
			for(i = 0; i < NUm_P; i += 1)
			pEx[][][i] = real(exp(sqrt(-1)*(2*pi/Num_P*i))*NormalEx[p][q][Num_wl])
			endfor
	Fullpath		= FolderPath+N_Ey+".mat"
	loadName = "N_"+N_Ey
		MLLoadWave /O/M=2/Y=4/S=2/N=$(LoadName) FullPath
		Wave/C NormalEy = $(loadName)
			re_f = real(NormalEy)
			im_f = imag(NormalEy)
			magnitude += re_f^2+im_f^2
			for(i = 0; i < Num_P; i += 1)
			pEy[][][i] = real(exp(sqrt(-1)*(2*pi/Num_P*i))*NormalEy[p][q][Num_wl])
			endfor
	Fullpath		= FolderPath+N_Ez+".mat"
	loadName = "N_"+N_Ez
		MLLoadWave /O/M=2/Y=4/S=2/N=$(LoadName) FullPath
		Wave/C NormalEz = $(loadName)
			re_f = real(NormalEz)
			im_f = imag(NormalEz)
			magnitude += re_f^2+im_f^2
			for(i = 0; i < Num_P; i += 1)
			pEz[][][i] = real(exp(sqrt(-1)*(2*pi/Num_P*i))*NormalEz[p][q][Num_wl])
			endfor
			
			magnitude	= sqrt(magnitude)
			pmag		= sqrt(pEx^2+pEy^2+pEz^2)

Redimension/n = (Num_x-1,Num_y-1,Num_F) magnitude,pEx,pEy,pEz,pmag	//画像表示用のセル数調整
killwaves NormalEx,NormalEy,NormalEz
		
end