本发明属于电子信息技术领域，具体为一种星载合成孔径雷达的极化系统参数设计方法。
背景技术：
全极化合成孔径雷达系统获取的图像能反映地物不同的散射机理和极化特性。利用全极化合成孔径雷达图像能开展很多应用，如极化目标分解，地表分类，军事目标探测，植被高度反演等。这些应用将有助于人们了解生态气候环境变化、资源分布，保卫领土安全等。
虽然极化技术优势显著，但极化雷达系统获取的图像不可避免地会受到系统扰动参数的污染。极化系统参数主要指极化串扰，通道不平衡和系统噪声。为了准确地使用极化信息，极化扰动参数必须得到校正。已有的研究多关注极化数据校正方法，并假设校正后的数据满足后续的应用需求[1][2]。只有少数研究关注具体的极化数据应用校正需求[3][4]。而对系统设计者来说，理解终端应用的系统参数设计需求是非常重要的。
针对极化数据的几个典型应用，即极化目标分解，地表分类，军事目标探测和植被高度反演开展极化系统参数设计需求分析非常重要。这些应用能帮助人们了解地球资源分布，气候环境变化和人类活动轨迹。我们必须给出系统极化参数设计需求，以确保设计的合成孔径雷达系统能按应用需求提供数据服务。为了提出极化系统参数设计需求，我们旨在建立极化系统参数与具体终端应用的评价指标之间的关联模型和误差传递关系。
对于极化目标分解和地表分类应用，许多极化目标分解参数可以用于敏感性分析，如取向角，cameron单位圆，cloude分解的熵参数h，极化旋转角参数α，去取向参数u和v[5][6][7][8]。由目标分解参数构成的参数空间可以用于分类目标和地表。这些目标分解方法在非螺旋体假设下是相互等价的，相应的目标分解参数也可以相互转化[9]。
军事目标探测算法有很多，它们的原理是测量两个散射矩阵的相似性。常用的算法有yang方法[10]和marino方法[11]。yang方法通过测量两个目标mueller矩阵的相似性来探测目标。marino方法通过测量两个目标的相干性系数来探测目标。
森林是生态系统的重要组成部分，探测森林的生物量有助于了解气候环境变化，而森林的生物量和森林垂直结构分布密切相关。为了反演植被高度，可以将极化信息引入干涉处理来分离不同散射机理的有效相位中心，该方法称为极化干涉技术[12][13]。欧空局已经计划在2020年开展biomass任务，通过极化干涉测量来绘制全球森林生物量[14]。中国也计划发射twin-l卫星来获取森林高度。
极化干涉技术反演植被高度是基于植被的电磁散射模型进行反演的，一个广泛使用的森林电磁散射模型为randomvolumeoverground(rvog)[15][16]。该模型将植被建模为一个两层结构，包含随机取向粒子层和地表层，而且可以用植被高度、衰减系数等少量垂直结构参数来描述。最终根据rvog模型，可以建立植被干涉相干性与植被参数的关系[17][18]。
利用已知的电磁散射模型，从极化干涉测量中就可以反演植被高度。测量值复干涉相干性是极化基的函数，利用相干优化技术可以获得最佳的高度检测[19]。相干优化过程就是选择能使相干性幅度最大的极化基的线性组合。使用相干优化技术可以获得最优的体散射机理和地散射机理，从而获得最佳的高度检测。
在通过数值实验建立极化系统参数与终端应用的评价指标之间的数值传递关系时，需要在仿真的全极化合成孔径雷达图像上人为地添加极化系统误差。对于极化合成孔径雷达图像的仿真，映射投影算法(mappingandprojectionalgorithm)可以仿真多种地表类型，如森林，农田，城区，道路，河流等[20]。双向解析射线追踪算法(bidirectionalanalyticraytracing)可以仿真多种人造目标，如飞机、坦克、舰船等[21]。polsarprosim森林仿真器可以仿真植被的极化干涉数据，该仿真器已经集成到欧空局polsarpro教育工具箱中[22]。
在对随机的植被高度反演误差进行统计建模时，需要使用概率分布模型。常用的统计分布模型有瑞利分布，伽马分布，k分布，逆高斯分布等[23][24][25][26]。对于从观测数据中检测未知的模型参数，可以使用的方法有矩估计，最大似然估计，对数累积量法等[27][28]。为了说明检测的模型参数的合理性，需要定量评估建模精度。常用的评价准则有k-s距离，它是通过比较观测数据与统计分布模型的累积分布函数之间的偏差来评估建模精度的[29]。如果偏差值为零，可以认为统计分布函数具有很好的建模精度。
技术实现要素：
本发明的目的在于提出一种星载合成孔径雷达的极化系统参数设计方法。该方法能为系统设计者提供星载合成孔径雷达的极化系统参数设计需求，确保设计的合成孔径雷达系统能按应用需求提供数据服务。
本发明从合成孔径雷达系统设计者的角度，建立了具体终端应用的评价指标与一种或几种极化系统参数与之间的关联模型，即理论模型公式或数值传递关系，进而可以得到相应终端应用的极化系统参数设计需求。本发明的技术方案具体介绍如下。
一种星载合成孔径雷达的极化系统参数设计方法，针对极化目标分解，地表分类，军事目标探测和植被高度反演这四种终端应用，通过建立终端应用评价指标和极化系统参数之间的关联模型进行极化系统参数设计，关联模型包括了理论模型公式或数值传递关系；所述极化系统参数包括极化串扰、通道不平衡和系统噪声中的一种或几种；其中：对于极化目标分解应用，建立评价指标和极化串扰之间的关联模型公式；对于地表分类应用，建立评价指标和极化串扰之间的数值关系；对于军事目标探测应用，建立评价指标和极化串扰之间的数值关系；对于植被高度反演应用，不仅建立评价指标与极化串扰、通道不平衡和系统噪声之间的关联模型公式，还得到评价指标和这三个极化系统参数之间的数值传递关系，理论模型公式和数值关系相互验证；具体设计方法如下：
终端应用为极化目标分解时，其评价指标为取向参数p和散射机理参数g，采用基于偏导数的灵敏度理论分析，建立极化串扰δ1和取向参数扰动量|δp|、散射机理参数扰动量|δg|之间的关联模型，该理论模型公式如下：
其中:δ1为极化串扰，p1、k为中间变量。
进而得到极化目标分解应用的极化串扰需求；
终端应用为地表分类时，其评价指标为极化分解熵参数h，极化旋转角参数α，极化去取向理论中反映地表介电特征的参数u，反映地表散射机制的参数v，以及地表分类误差；四个极化分解参数h、α、u和v和地表分类误差对极化串扰δ2的敏感性采用映射投影算法mpa进行仿真包含多种地表类型的合成孔径雷达图像来进行敏感性分析，建立极化串扰δ2与h，α，u，v，以及地表分类误差之间的数值传递关系，进而得到地表分类应用的极化串扰需求；
终端应用为军事目标探测时，其评价指标为目标探测精度，目标探测精度对极化串扰δ3的敏感性采用双向解析射线追踪算法bart进行仿真多种场景来进行敏感性分析，建立极化串扰δ3与目标探测精度之间的数值传递关系，进而得到军事目标检测应用的极化串扰需求；
终端应用为植被高度反演时，其评价指标为植被高度反演误差δ，对由极化系统参数构成的扰动矩阵进行特征值分析，建立相干优化法植被高度反演误差δ与极化串扰δ4、通道幅度不平衡|f|、噪声功率
其中：δ4为极化串扰，
此外，终端应用为植被高度反演时，随机植被高度反演误差的均值对极化串扰δ4、通道幅度不平衡|f|和噪声功率
本发明中，所述终端应用为地表分类时，极化参数组成的u-v-h参数空间对不同的地表类型进行分类；地表在参数空间u-v-h上被分为9类，通过统计错误分类的像素得到分类误差。
本发明中，所述终端应用为军事目标探测时，目标探测方法采用yang方法和marino方法。
本发明中，所述终端应用为植被高度反演时，对仿真的森林场景人为添加不同的极化系统误差，进而对随机的植被高度检测误差用gamma分布进行建模，检测未知的模型参数μ可以得到反演误差的均值。
和现有技术相比，本发明的有益效果在于：本发明从系统设计者的角度出发，针对极化数据的几个典型应用开展了极化系统参数设计需求分析。这些应用能帮助人们了解地球资源分布，气候环境变化和人类活动轨迹。我们给出的系统极化参数设计需求，可以确保设计的合成孔径雷达系统能按应用需求提供数据服务。本发明为星载合成孔径雷达的系统设计者提供了极化系统参数设计方法。
附图说明
图1为本发明的实施例中散射机理扰动量|δg|在cameron单位圆上的幅度图。
图2为本发明的实施例中取向参数扰动量|δp|和散射机理扰动量|δg|在g-p平面上的幅度变化图。
图3为本发明的实施例中基于mpa算法仿真得到的合成孔径雷达图像。
图4为本发明的实施例中不同地表类型区域在u-v-h平面上的分布。
图5为本发明的实施例中极化分解参数h、α、u和v随极化串扰δ2的变化曲线。
图6为本发明的实施例中地表分类误差随极化串扰δ2的变化曲线。
图7为本发明的实施例中sandiego港口的uav合成孔径雷达图像和谷歌地球航拍图像。
图8为本发明的实施例中sandiego港口的极化分解参数h、α、u和v随极化串扰δ2的变化曲线。
图9为本发明的实施例中三种典型军事目标模型。
图10为本发明实施例中基于bart算法仿真得到的不同军事目标的合成孔径雷达图像。
图11为本发明实施例中yang探测器的探测误差随极化串扰δ3的变化曲线。
图12为本发明实施例中marino探测器的探测误差随极化串扰δ3的变化曲线。
图13为本发明的实施例中polsarprosim森林场景仿真示意图和pauli矢量伪彩色编码合成孔径雷达图像。
图14为本发明的实施例中的森林高度反演结果和森林高度分布直方图。
图15为本发明的实施例中gamma分布拟合植被高度反演误差的概率分布直方图和累积分布直方图的结果。
图16为本发明的实施例中gamma分布函数建模精度评估。
图17为本发明的实施例中高度检测误差理论公式和参数μ随极化串扰δ4的变化结果对比。
图18为本发明的实施例中高度检测误差理论公式和参数μ随幅度不平衡|f|的变化结果对比。
具体实施方式
下面结合附图和实施例对本发明的技术方案进行详细介绍。
本发明提出的星载合成孔径雷达的极化系统参数设计方法，从合成孔径雷达系统设计者的角度，建立了具体终端应用的评价指标与极化系统参数与之间的关联模型。该关联模型包括了理论模型公式和数值传递关系。基于该关联模型可以提出相应终端应用的极化系统参数设计需求。下面结合具体的实施例对本发明做详细说明。
对极化目标分解应用，本发明建立了极化串扰δ1与极化目标分解的取向参数p和散射机理参数g之间的关联模型理论公式。极化分解参数可以表示为极化串扰的函数p(δ1)和g(δ1)。将函数p(δ1)和g(δ1)对δ1求一阶偏导数可以得到极化分解参数扰动量|δg|和|δp|与极化串扰δ1之间的关联模型公式。该关联模型适用于确定性的相干散射目标。这两个极化分解参数区分了极化散射目标的不同类型，其数值可以在cameron单位圆上展示。
在实施例中，我们根据推导得到的关联模型公式研究了扰动量|δg|和|δp|的模式图。散射机理扰动量|δg|在cameron单位圆上的梯度图如图1所示。图1表明散射机理参数g的主要偏移是从边缘到中心，这种收缩偏移和取向角p以及散射机理参数g是成比例的。除了主要偏移，当取向角很小时还可以观察到一种次级偏移在区域g→1时产生。此外，扰动量|δg|和|δp|在参数g-p平面上的偏移如图2所示。此时，散射机理参数g的收缩偏移仍然能清晰观察到。取向的偏移主要发生在右下角区域，它会使得该区域目标的取向角增大。同时，左下角区域不受极化串扰的影响，也就是极化特性对串扰不敏感。
对地表分类应用，本发明建立了极化串扰δ2与极化分解熵参数h，极化旋转角参数α，地表介电参数u，地表散射机制参数v，以及地表分类误差之间的数值传递关系。在实施例中，我们使用映射投影算法mpa仿真了l波段分辨率为12m的合成孔径雷达图像，仿真场景包含森林，农田，城区，道路，河流等多种地表类型，仿真得到的图像如图3所示。在图3中，我们用红框选择了多种地表类型用于后续的敏感性分析。本实施例首先计算了选择区域的极化分解参数，并在参数空间u-v-h上对这些地表进行分类，分类结果如图4所示。从图4看出，这些区域在u-v-h参数空间上得到了有效分离。
对仿真的合成孔径雷达图像添加不同的极化串扰并对所选区域进行极化分解，分析四个极化分解参数h、α、u和v对极化串扰δ2的敏感性。其中，极化串扰δ2从-42db变化到-18db。在本实施例中，四个极化分解参数随极化串扰的变化曲线如图5所示。从图中可以看出，当极化串扰小于-32db时，四个极化分解参数的变化可以忽略。参数u的变化比其它三个参数更敏感，当极化串扰小于-25db时，参数u的最大偏移小于0.1，所以我们认为对于自然地表分类，极化串扰满足小于-25db的要求即可。
本实施例中，在计算得到典型地表区域的极化分解参数后，可以在参数空间u-v-h上将地表分为9类。在不同的极化串扰下，在参数平面u-v-h上统计错误分类的像素可以得到分类误差。分类误差随极化串扰δ2的变化曲线如图6所示。从图中看出，当极化串扰小于-25db时，分类误差低于10％。而当极化串扰从-25db增大到-18db时，分类误差迅速从10％增大到44％。这一结果也证实了我们对自然地表分类提出的-25db极化串扰需求。通过本实施例得到的地表分类应用的极化串扰设计需求为：为使得极化分解参数不随极化串扰值而改变，且自然地表分类误差小于10％，应保证极化串扰小于-25db。
对地表分类应用，本发明还在实施例中对nasa/jpl获取的sandiego港口的uav合成孔径雷达图像进行了极化串扰敏感性分析。uav合成孔径雷达图像以及对应的谷歌航拍图像如图7所示。我们在图中标记了三块区域进行分析，分别是城区、海面和森林。我们假设该数据不受极化串扰污染，通过仿真重新添加极化串扰。在本实施例中我们分析了三块区域的四个极化分解参数h、α、u和v随极化串扰δ2的变化曲线，结果如图8所示。从图8中我们可以得到同样的结论，即自然地表分类应满足至少-25db的极化串扰需求。
对军事目标探测应用，本发明建立了极化串扰δ3与目标探测误差之间的数值传递关系。在实施例中，我们使用双向解析射线追踪算法bart仿真了多种场景来进行敏感性分析，仿真的典型场景包括机场上的飞机，草地上的坦克，以及海面上的舰船。本实施例中仿真的飞机、坦克、舰船模型如图9a-9c所示。我们仿真了机场上的15架飞机，草地上的30架坦克和海面上的30艘舰船，这些目标共具有0°,45°,90°,135°度四种方位角。仿真得到的合成孔径雷达图像具有1m分辨率，仿真结果如图10a-10c所示。
在本实施例中，我们使用yang方法和marino方法进行目标探测，两者都是通过测量散射矩阵的相似性来探测目标的。yang方法测量mueller矩阵相似性，marino方法测量散射矢量的相干性系数。对仿真场景添加不同的极化串扰，使用两种探测方法分别进行目标探测，可以得到目标探测误差随极化串扰δ3的变化曲线。yang探测器的探测误差随极化串扰δ3的变化曲线如图11a-11c所示，marino探测器的探测误差随极化串扰δ3的变化曲线如图12a-12c所示。从变化曲线可以看出，当极化串扰大于-25db时，探测误差很快增大。通过本实施例得到的军事目标探测应用的极化串扰设计需求为：对于仿真的典型场景，为获得至少90％的目标探测精度，应保证极化串扰小于-25db。
对植被高度反演应用，本发明建立了植被高度反演误差δ与极化串扰δ4、幅度不平衡|f|以及噪声功率
在本实施例中，对仿真得到的原始森林场景使用相干优化技术反演植被高度，反演结果将作为森林的高度真值。对每个像素求解最优相干性并检测高度，检测的森林高度分布图和森林高度分布直方图如图14所示。森林区域在图中凸显出来，且高度基本一致。从高度分布直方图中看出，森林高度以10m为中心，分布在较小范围内。统计森林高度的均值为10.06m，和仿真设置的参数很接近，说明反演方法效果很好。
在本实施例中，对仿真得到的森林场景人为添加不同的极化系统误差，并对随机的植被高度反演误差用gamma分布进行建模。以极化串扰δ4＝0.1，snr＝20db为例，从受污染的数据中检测森林高度并统计每个像素的反演误差。使用gamma分布拟合反演误差并检测未知的模型参数。图15给出了gamma分布拟合概率分布直方图和累积分布直方图的结果。图中蓝条代表观测到的检测误差，红色拟合曲线为具有特定模型参数的gamma分布函数。从图15中可以看到，gamma分布函数能很好地匹配观测数据，也就是能给植被反演误差提供很好的建模能力。
在本实施例中，gamma分布具有两个模型参数μ和l，参数μ代表分布的均值，参数l影响gamma分布的形状。通过调节参数l可以获得最佳的拟合曲线，参数μ反映了数据的特征，所以我们统计了参数μ与极化系统参数之间的数值传递关系。对不同信噪比和极化串扰组合下的高度检测误差进行gamma分布建模并检测一组模型参数μ。同理，对不同信噪比和幅度不平衡组合下的高度检测误差进行gamma分布建模并检测一组模型参数μ。
在本实施例中，为了说明检测的模型参数是合理的，我们使用k-s距离来评估建模精度。对不同极化系统参数组合下的建模精度用k-s距离来评估，结果如图16所示。从图中看出，k-s距离都非常小，这说明我们使用的gamma分布模型对观测数据的拟合效果好，偏差小。好的建模精度说明我们在不同极化系统参数组合下检测的模型参数μ是正确的。
在本实施例中，我们得到了模型参数μ和极化串扰、通道幅度不平衡和信噪比之间的数值传递关系。该数值关系与本发明中建立的植被高度反演误差δ与三个极化系统参数之间的关联模型理论公式是一致的，两者可以相互验证。在图17中，我们对比了高度检测误差理论公式和参数μ随极化串扰δ4的变化。在图18中，我们对比了高度检测误差理论公式和参数μ随幅度不平衡|f|的变化。在两幅图中，离散点代表从观测数据中检测到的参数μ值，虚线代表本发明中建立的关联模型公式。从图中可以看出，虚线和离散点基本吻合，且偏差非常小。这说明我们建立的植被高度检测误差和极化系统参数之间的关联模型理论公式和数值传递关系是一致的。通过本实施例得到的植被高度反演应用的极化系统参数设计需求为：为使得由极化系统误差引起的平均植被高度反演误差小于0.5m，信噪比应当优于18db，极化串扰应当小于-25db，幅度不平衡应当小于0.5db。
参考文献
[1]j.j.vanzyl,calibrationofpolarimetricradarimagesusingonlyimageparametersandtrihedralcornerreflectorresponses，ieeetrans.geosci.remotesens.,vol.28,no.3,pp.337-348,may.1990.
[2]a.freeman,j.j.vanzyl,j.d.klein,h.a.zebker,andy.shen，calibrationofstokesandscatteringmatrixformatpolarimetricsardata，ieeetrans.geosci.remotesens.,vol.30,no.3,pp.531-539,may.1992.
[3]y.wangandv.chandrasekar,polarizationisolationrequirementsforlineardual-polarizationweatherradarinsimultaneoustransmissionmodeofoperation,ieeetrans.geosci.remotesens.,vol.44,no.8,pp.2019-2028,2006.11.
[4]r.touzi,p.w.vachon,andj.wolfe,requirementonantennacross-polarizationisolationfortheoperationaluseofc-bandsarconstellationsinmaritimesurveillance,ieeegeosci.remotesens.lett.,vol.7,pp.861-865,oct.2010.
[5]w.l.cameronandl.k.leung,featuremotivatedpolarizationscatteringmatrixdecomposition,ieee1990int.radarconf.,pp.549-557,1990.
[6]j.r.huynen,phenomenologicaltheoryofradartargets.phd.thesis,univ.technologydelft,netherlands,1970.
[7]s.r.cloudeande.pottier,anentropybasedclassificationschemeforlandapplicationofpolarimetricsar,ieeetrans.geosci.remotesens.,vol.35,no.1,pp.68-78,1997.
[8]f.xu；y.-q.jin,deorientationtheoryofpolarimetricscatteringtargetsandapplicationtoterrainsurfaceclassification,ieeetrans.geosci.remotesens.,vol.43,no.10,pp.2351,2364,2005.
[9]y.q.jin,andf.xu,polarimetricscatteringandsarinformationretrieval,wiley-ieee,2013.
[10]j.yang,y.n.peng,ands.m.lin,similaritybetweentwoscatteringmatrices.electron.lett.,vol.37,no.3,pp.193-194,2001.
[11]a.marino,s.r.cloude,andi.h.woodhouse,detectingdepolarizedtargetsusinganewgeometricalperturbationfilter,ieeetransactionsongeoscienceandremotesensing,vol.50,no.10,pp.3787-3799,2012.
[12]s.r.cloudeandk.p.papathanassiou,polarimetricsarinterferometry,ieeetrans.geosci.remotesens.,vol.36,no.5,pp.1551-1565,sep.1998.
[13]k.p.papathanassiouands.r.cloude,single-baselinepolarimetricsarinterferometry,ieeetrans.geosci.remotesens.,vol.39,no.11,pp.2352-2363,nov.2001.
[14]t.l.toan,ap-bandsarforglobalforestbiomassmeasurement:thebiomassmission,ingeneralassemblyandscientificsymposium(ursigass),xxxithursi,aug.2014.
[15]r.n.treuhaft,s.n.madsen,m.moghaddam,andj.j.vanzyl,"vegetationcharacteristicsandunderlyingtopographyfrominterferometricradar,"radiosci.,vol.31,no.6,pp.1449-1495,nov.1996.
[16]r.n.treuhaftandp.r.siqueira,"verticalstructureofvegetatedlandsurfacesfrominterferometricandpolarimetricradar,"radiosci.,vol.35,no.1,pp.141-177,jan.2000.
[17]i.hajnsek,f.kugler,s.k.lee,andk.p.papathanassiou,"tropicalforestparameterestimationbymeansofpol-insar:theindrex-iicampaign,"ieeetrans.geosci.remotesens.,vol.47,no.2,pp.481-493,feb.2009.
[18]s.r.cloudeandk.p.papathanassiou,"three-stageinversionprocessforpolarimetricsarinterferometry,"proc.inst.elect.eng.—radarsonar,navig.,vol.150,no.3,pp.125-134,jun.2003.
[19]s.r.cloudeandk.p.papathanassiou,"polarimetricsarinterferometry,"ieeetrans.geosci.remotesens.,vol.36,no.5,pp.1551-1565,sep.1998.
[20]f.xuandy.q.jin,"imagingsimulationofpolarimetricsarforacomprehensiveterrainsceneusingthemappingandprojectionalgorithm,"ieeetrans.geosci.remotesens.,vol.44,no.11pp.3219-3234,2006.
[21]f.xuandy.q.jin,bidirectionalanalyticraytracingforfastcomputationofcompositescatteringfromanelectric-largetargetoverrandomlyroughsurface,ieeetransactionsonantennasandpropagation,2009,57(5):1-11.
[22]m.l.williamsande.pottier,forestcoherentsarsimulationwithinpolsarpro:aneducationaltoolboxforpolsarandpolinsardataprocessing,inasianconferenceonremotesensing,kualalumpur,malaysia,nov,2007.
[23]j.w.goodman,statisticalpropertiesoflaserspecklepatterns:springerberlinheidelberg,1975.
[24]c.oliverands.quegan,understandingsyntheticapertureradarimages:scitechpublishing,2004.
[25]e.jakemanandp.n.pusey,amodelfornon-rayleighseaecho,ieeetrans.antennaspropagat.,vol.24,no.6,pp.806-814,nov.1976.
[26]a.c.frery,h.j.muller,c.c.f.yanasse,ands.j.s.sant'anna,amodelforextremelyheterogeneousclutter,ieeetrans.geosci.remotesensing,vol.35,no.3,pp.648-659,may,1997.
[27]m.g.kendall,a.stuart,andj.k.ord,kendall'sadvancedtheoryofstatistics:universitypress,1994.
[28]v.a.krylov,g.moser,s.b.serpico,andj.zerubia,"onthemethodoflogarithmiccumulantsforparametricprobabilitydensityfunctionestimation,"ieeetrans.imageprocess.,vol.22,no.10,pp.3791-3806,oct.2013.
[29]m.d.devoreandj.a.o'sullivan,"quantitativestatisticalassessmentofconditionalmodelsforsyntheticapertureradar,"ieeetrans.imageprocess.,vol.13,no.2,pp.113-125,feb.2004.