第3章

数据分析处理



本章包含22个案例，介绍利用Igor进行数据分析和处理的程序设计技巧。内容涵盖： ①数据分析处理中wave操作的程序设计； ②特定领域里典型数据处理分析的程序设计。
数据分析处理程序设计一般包括数据拟合、算法设计、数据管理、数据可视化等操作，个别情况下还需要考虑程序运行效率。程序形式既有命令行程序，也有图形界面程序，具有较强的综合性，对程序设计者的要求是比较高的。
一般而言，好的数据分析处理程序应有以下特点： ①函数功能明确，一个函数只解决一个问题； ②同一个功能尽可能封装成为一个函数，在程序里需要该功能时调用该函数，而不是再将所有的代码复制一遍； ③尽可能不要在完成运算功能的函数里进行图形显示或者创建数据文件夹等操作。
案例43对wave的基本操作
1. 应用情景
对wave的基本操作包括获取wave名字、赋值、设置坐标、获取坐标信息、显示、显示设置等。
2. 程序代码

//******************************************************

//1. 对wave的基本操作

//******************************************************

Function WaveOPDaemon()



//创建wave

Make/O ddd



//设置wave坐标

SetScale/P x,0,1,"s",ddd//设置横坐标

SetScale d,0,0,"m",ddd    //设置数据单位

SetDimLabel 0,1,Row1,ddd//设置行名



//赋值

ddd=sin(x)



//操作

Printf "sum:%f＼r", sum(ddd)  //求所有数据和

Printf "area:%f＼r" area(ddd)   //求曲线包围面积

Printf "average:%f＼r" sum(ddd)/numpnts(ddd) //求平均值

Printf "ddd［20］:%f＼r" ddd［20］	    //显示第21个数据(从0开始计数)

Printf "ddd［Row1］:%f＼r" ddd［%Row1］ //显示行名为Row1所在行数据

Duplicate/O ddd ddd1         //复制数据

ddd1=ddd+cos(x)              //对ddd1利用ddd进行赋值

Edit ddd	                    //用表格显示数据

Edit ddd.xy	                 //用表格显示数据及其x坐标值

Edit ddd.x                   //用表格显示x坐标值

Edit ddd.ld						 //用表格显示数据和行名

Edit ddd.l						 //用表格显示行名

Display ddd,ddd1             //显示曲线

ModifyGraph mode(ddd)=3,marker(ddd1)=8



//获取数据信息

Printf "Size:%d＼r" numpnts(ddd)　　//获取数据点数

Printf "x0:%f＼r" leftx(ddd)   //获取第一个数据点x坐标值

Printf "x0:%f＼r" DimOffset(ddd,0) //获取第一个数据点x坐标值

Printf "Size:%d＼r" DimSize(ddd,0)　　//获取数据点数

Printf "delta x:%f＼r" DimDelta(ddd,0) //获取x坐标间隔

Print NameofWave(ddd) //获取wave名字

End





3. 程序分析
Make命令创建wave，SetScale命令设置wave的坐标，标记/I指定wave首末数据点对应的坐标值，其他点的坐标值自动计算(均匀分布)，也可使用标记/P，对应的命令如下： 

SetScale/P x,x0,dx,ddd

其中： x0表示第一个数据点对应的坐标值； dx表示坐标间隔。
SetScale d,0,0,"m",ddd中的d表示设定数据的单位为m。
SetDimLabel 0,1,Row1,ddd设置行名为Row1。0表示设置行名，1表示设置行名的行是第1行。SetDimLabel也可以设置列名。关于SetDimLabel的详细使用请查看帮助文件。
ddd=sin(x)利用x函数对wave进行赋值。整个赋值过程是自动的： 先取出某个数据点，计算该点的x坐标值，利用sin函数计算该坐标值的正弦值，将值赋值给该数据点，再循环对下个数据点进行赋值。因此这条命令在Igor内部是一个循环过程。上面的赋值过程可以利用循环改写如下： 

Variable i,n=numpnts(ddd)

Variable x0=leftx(ddd)

Variable dx=deltax(ddd)

Variable xx

for(i=0;i<n;i+=1)

xx=x0+dx*i

ddd［i］=sin(xx)

endfor

可见利用ddd=sin(x)可大大简化赋值过程。
sum函数计算所有数据点的和，area函数计算曲线和x坐标轴所包围的面积。注意这两个函数是不一样的，前者是将所有数据点加起来，不考虑x坐标值； 后者相当于计算曲线积分值，考虑了x坐标，area(ddd)=sum(ddd)*deltax(ddd)。
ddd［20］利用中括弧获取第21个点的数据。ddd(0.1)利用小括弧获取坐标值为0.1时该点的数值，如果这个坐标不存在，则为插值结果，这是最简单的一种插值方法。利用小括弧也可以设置某个坐标处的数值，如ddd(0.1)=sin(0.1)。
ddd［%Row1］获取行名为Row1的行的值。



注意： 
行名和列名还能够描述行或者列的功能，如绘制风羽图时，箭头信息wave的第3列列名必须指定为WindBarb，Igor才能够将该列正确识别为风速大小列； 设置自定义刻度标签时，存放刻度标签的wave第2列列名必须指定为Tick Type，Igor才能够将该列识别为刻度线列。


Duplicate复制一个wave，复制的内容包括数据和坐标信息，因此是完全的复制。在复制时可以指定范围，如只复制坐标值区间(0.1,0.2)内的内容，则命令为： 

Duplicate/O/R=(0.1,0.2) ddd ddd1

如果将小括弧换成中括弧，则表示指定复制的数据点范围。
ddd1=ddd+cos(x)给ddd1赋值，赋值过程仍然是一个循环。
Edit ddd利用表格显示wave。此命令只显示wave的数据，如果要显示x坐标，需要使用Edit ddd.xy。
Edit ddd.ld和Edit ddd.l分别用表格显示行名/数据或行名。
Display ddd显示wave。
numpnts计算wave的数据点数，leftx获取wave第一个数据点x坐标值，delta获取坐标间隔。这3个函数要求wave的维度是一维的。
DimSize和numpnts含义一样，DimOffset和leftx含义一样，DimDelta和deltax含义一样，只是前者适用范围更广，可以获取任意维数wave的坐标信息。调用这3个函数时需要指定wave的维数，0表示行(x坐标)，1表示列(y坐标)，2表示层(z坐标)。
NameOfWave得到wave的名字。 Igor程序设计中很多命令需要使用wave的名字作为参数，通过此函数就可以获取wave名字。
4. 实战应用
这是一个功能演示程序，在命令行窗口直接输入WaveOPDaemon()，观察执行结果，如图31所示。


图31对wave的基本操作


案例44获取wave的坐标信息
1. 应用情景
wave是Igor最核心的概念，坐标是wave的重要属性。数据处理中经常需要获知wave的坐标信息，坐标信息包括起始值、终止值、坐标间隔以及数据点数等。下面的函数封装了获取wave坐标信息的操作，程序只需要调用此函数就可以获取wave的坐标信息。
2. 程序代码

//*************************************************************

//1.获取wave的坐标信息

//************************************************************* 

Function GetDimInfo(w,nx,dx,x0,ny,dy,y0)

Wave w//2d image

Variable &nx//行

Variable &dx

Variable &x0

Variable &ny//列

Variable &dy

Variable &y0

nx=DimSize(w,0)

dx=DimDelta(w,0)

x0=DimOffset(w,0)

ny=DimSize(w,1)

dy=DimDelta(w,1)

y0=DimOffset(w,1)

End

3. 程序分析
程序中调用坐标属性函数获取wave的坐标信息，然后赋值给相应的变量，nx表示行数(x坐标方向数据点数)，dx表示行坐标间隔(x坐标间隔)，x0表示最小的行坐标(x坐标值)，ny表示列数(y坐标方向数据点数)，dy表示列坐标间隔(y坐标间隔)，y0表示最小的列坐标(y坐标值)。使用此程序可以很方便地获取wave的坐标属性，无须多次重复调用DimOffset等函数。



小技巧： 
程序里使用了变量引用作为函数参数的技巧，通过传递变量的引用，可以使函数返回多个值。


4. 实战应用
在程序窗口输入以下程序： 

Function Daemon()

Variable nx,dx,x0,ny,dy,y0

Make/O/N=(20,20) ddd

SetScale/P x,0,0.1,ddd

SetScale/P y,0,0.3,ddd

GetDiminfo(ddd,nx,dx,x0,ny,dy,y0)

Print nx,dx,x0,ny,dy,y0

End

程序创建一个名为ddd的wave，利用GetDimInfo获取ddd的坐标信息，并输出到命令行窗口，结果如图32所示。


图32获取wave的坐标信息

案例45对二维wave按照行或者列归一化
1. 应用情景
由于仪器测量通道的不均匀性，数据中会存在条纹状的结构，这可以通过按照行或列对数据进行归一化消除。
2. 程序代码

//*************************************************************

//1. 对二维wave按照行或者列归一化

//************************************************************* 

Function NormalizeImage(w,dim,v1,v2)

Wave w　　//Image

Variable dim　　//归一化方向,1:x,2:y

Variable v1,v2 //归一化起始和终止坐标值

Variable i,j,nx,ny,dx,dy,x0,y0,sum0

GetDimInfo(w,nx,dx,x0,ny,dy,y0)

if(1==dim)

ny=DimSize(w,1)

for(j=0;j<ny;j+=1)

MatrixOP/O out=col(w,j)

SetScale/P x,x0,dx,out

sum0=sum(out,v1,v2)

w［］［j］=w［p］［j］/sum0

endfor

endif

if(2==dim)

nx=DimSize(w,0)

for(i=0;i<nx;i+=1)

MatrixOP/O out=row(w,i)

Redimension/N=(numpnts(out)) out

SetScale/P x,y0,dy,out

sum0=sum(out,v1,v2)

w［i］［］=w［i］［q］/sum0

endfor

endif	

End

//*************************************************************

//2. 获取wave坐标信息

//************************************************************* 

Function GetDiminfo(w,nx,dx,x0,ny,dy,y0)

Wave w//2d image

Variable &nx　　//行

Variable &dx

Variable &x0

Variable &ny　　//列

Variable &dy

Variable &y0

nx=DimSize(w,0)

dx=DimDelta(w,0)

x0=DimOffset(w,0)

ny=DimSize(w,1)

dy=DimDelta(w,1)

y0=DimOffset(w,1)

End

3. 程序分析
程序由两部分组成： 第一部分完成归一化功能； 第二部分作为辅助函数获取wave的坐标信息。
第一部分函数NormalizeImage实现归一化功能，有4个输入参数： w表示待归一化wave； dim取值为1对行归一化，取值为2对列归一化； v1和v2设定获取归一化因子的范围。归一化的基本原理是： 按行或者列抽出每一条曲线，计算坐标值在v1和v2之间所有的数据点的和sum0，然后利用sum0除整条曲线。程序首先利用MatrixOP命令抽取一行或者一列，然后利用SetScale命令设置曲线的坐标，再利用sum函数计算坐标值落在v1和v2之间的所有数据点的值的和，最后利用wave的自动赋值技巧将该行或者列所有的值除以和值。抽取行或者列也可以通过别的方法进行，但是利用MatrixOP命令是最方便、最有效的。
第二部分函数GetDimInfo介绍见本章案例44。
4. 实战应用
假设待归一化的wave名为AuLaser，位于当前目录下，归一化前如图33所示。


图33直接测量的二维谱（存在通道不均匀性）


在命令行窗口输入： 

NnormalizeImage(AuLaser,1,2.56,2.6)

归一化后结果如图34所示。


图34按行(通道)归一化后的二维谱


案例46按行或者列获取二维谱所有曲线
1. 应用情景
对ARPES谱进行EDC或者MDC分析时，需要从二维谱中抽取某个能量值对应MDC或者某个动量值对应的EDC。在进行色散拟合时需要获取构成二维谱的所有MDC。为了改善数据统计，有时需要将相邻的几条曲线数据加在一起，以改善统计。下面的程序演示了如何从二维谱中抽取曲线。程序能够根据条件自动创建所有曲线，并返回曲线的名字列表，以供接下来使用。
2. 程序代码

//************************************************************************

//按照行或者列抽取二维谱所有曲线

//************************************************************************

Function/s GetTraces(w,nb,dim,［p0,p1,dp,x1,x2,df,xsym,symm factor］)

Wave w	　　//数据wave

Variable nb　　//Binding 数

Variable p0,p1,dp//起始点，终止点，间隔(同时提供)

Variable x1,x2   //坐标间隔

Variable dim     //0:沿行抽取 1:沿列抽取

String df        //抽取曲线存放于df中

Variable xsym,symm,factor

String curr=GetDataFolder(1)

if(numtype(strlen(df))!=2)

SetDataFolder df

endif

Variable i,j,nx,ny,dx,dy,x0,y0,k=0

String s,s1="",sw

if(dim==1)

Matrixtranspose w

endif

sw=NmeOfWave(w)

if(!ParamIsDefault(x1))

Duplicate/O/R=(x1,x2) w,tmp

Wave w=$nameofwave(tmp)

endif

GetDimInfo(w,nx,dx,x0,ny,dy,y0)

ny=dimsize(w,1)

if(ParamIsDefault (p0))

p0=0

p1=ny

dp=1

endif

if(ParamIsDefault (symm))

symm=0

endif

p1=p0+round((p1-p0)/dp)*dp

for(i=p0;i<=p1;i+=dp)

k=0

s=sw+"_"+num2str(i)

s1=s1+s+";"

Make/O/N=(nx) $s=0

SetScale/P x,x0,dx,$s

Wave w1=$s

for(j=round(-nb/2);j<=round(nb/2);j+=1)

if((i+j)<0||(i+j)>=ny)

continue

endif

MatrixOP/O out=col(w,i+j)

w1=w1+out

k+=1

endfor

w1=w1/k

endfor

if(dim==1)

Matrixtranspose w

endif

if(!ParamIsDefault (x1))

KillWaves/Z tmp,out

endif

SetDataFolder curr

s1=s1［0,strlen(s1)-2］

return s1

End

3. 程序分析
本程序参数较多，有3个输入参数，9个可选参数。
参数w表示要抽取曲线的二维wave； nb表示积分参数(即相邻多少条曲线数据加在一起)； dim表示抽取曲线的方向，dim取0时沿x方向抽取曲线，取1时沿y方向抽取曲线。
可选参数p0、dp和p1表示抽取曲线的范围，即从第p0条开始每隔dp条抽取一条曲线，一直到p1为止； x1和x2的含义一样，只是表示坐标而不是曲线次序； 如果不提供这些参数则抽取范围默认为整个wave，间隔为1或者最小坐标间隔。
df是文件夹名字，表示新创建的曲线将存放在此文件夹内，如果不提供则存放在当前文件夹，factor是采样因子。
程序虽然较长，但是原理并不复杂，基本思路是利用MatrixOP命令循环取出每一条曲线，再利用一个循环将相邻的nb条曲线和该曲线数据加在一起。当dim等于1时程序先利用MatrixTranspose命令交换w的行和列，然后再执行下面的操作。
程序里调用了案例44 GetDimInfo函数获取w的坐标信息。
4. 实战应用
假设被操作的二维wave为Bi2212，且位于当前目录下。打开程序窗口，输入下面的程序： 

Function Exp21()

Wave w=Bi2212

String s0=GetTraces(w,10,0)

Print s0

Wave EDC1=$StringFromList(0,s0)

Display EDC1

End

程序按照默认方式调用GetTraces函数，从第0条EDC开始，抽取所有EDC，每条EDC都是相邻10条EDC叠加在一起的结果。生成的所有EDC位于当前数据文件夹目录下。程序打印了每条EDC的名字，并显示第一条EDC。程序执行结果如图35所示。