MATLAB之TOPAS輸出檔案擬合高斯分佈並提取引數
阿新 • • 發佈:2018-12-31
此處需要將TOPAS跑出來的資料進行高斯擬合,並計算標準差sigma,以準備下一步的分析使用。csv資料型別如上篇文章所述,因為對TOPAS程式直接設定了輸出檔案的格式是“能量MeV+束斑大小cm+RangeShift寬度x空氣間隙寬度”,此處用迴圈讀取程式並輸出到一個儲存結果的xlsx檔案中。本來可以直接用list=dir()函式以及list(i).name來跑程式,可是給我的檔案命名格式有點問題,所以得重寫迴圈。MATLAB用dir讀取檔名的時候給出的排序都是從第一位開始依次比較大小,所以一開始定好檔名格式對後面的程式幫助很大。
clc; clear; fileDir='C:\Users\32628\Documents\MATLAB\results\'; %建立用於存放Sigma的陣列 valueM=zeros(10,18); row=1; column=1; for k=70:10:100 for j=1:1:3+(k-70)/10 for i=5:5:50 %讀取檔名 kk=num2str(k); jj=num2str(j); ii=num2str(i); filename=strcat(kk,'MeV',jj,'cm',ii,'.csv'); filename2=strcat(fileDir,filename); file=csvread(filename2,8,0); yvalue=file(:,4); xvalue=file(:,1); %設定要用的擬合模型 [xData, yData] = prepareCurveData( xvalue, yvalue ); ft = fittype( 'gauss1' ); opts = fitoptions( 'Method', 'NonlinearLeastSquares' ); %opts.Display = 'Off'; %opts.Lower = [-Inf -Inf 0]; %opts.StartPoint = [39194592 299 24.6030834338969]; %將模型用於資料 [fitresult, gof] = fit( xData, yData, ft, opts ); % Plot fit with data. %figure( 'Name', 'Normal Fit' ); %h = plot( fitresult, xData, yData ); %legend( h, 'yvalue vs. xvalue', 'untitled fit 1', 'Location', 'NorthEast' ); % Label axes %xlabel xvalue %ylabel yvalue %grid on valueM(row,column)=(fitresult.c1)/2^(1/2); if i~=50 row=row+1; else row=1; column=column+1; end end end end xlswrite('C:\Users\32628\Documents\MATLAB\sigma.xlsx',valueM,'Sheet1','B2');
本來想直接用normfit()來得到sigma的,但發現用這個函式得到的sigma有問題,才另想辦法用cftool來擬合數據並用Generate Code來跑。