资源描述:
《InSAR干涉纹图噪声抑制方法研究》由会员上传分享,免费在线阅读,更多相关内容在学术论文-天天文库。
InSAR干涉纹图噪声抑制方法研究摘要干涉合成孔径雷达(InSAR)三维成像技术具有覆盖面积大、空间分辨率高、高程精度高的优点,并且可以全天时、全天候地工作,是获取三维数字地面高程模型((DEM)的一种有效方法,可以广泛应用于地形测绘、地震和火山监测以及地形沉降检测等领域。论文在总结了InSAR三维成像技术的发展历史和现状后,分析了InSAR系统的工作原理和影响InSAR系统测量精度的各种因素,给出了系统的数据处理流程。InSAR技术是根据两幅SAR复图像的干涉相位差来获取地面目标的三维信息,因此干涉纹图中的噪声是影响系统测量精度的重要因素之一。为了获得准确的干涉相位差数据,使二维相位解缠能够顺利进行,必须对各种噪声进行有效的抑制,所以有必要对InSAR干涉相位纹图的噪声抑制方法进行研究。论文在研究InSAR干涉纹图噪声来源、特性和相位的统计特征的基础上,基于实验结果的基础上分析了传统的干涉纹图滤波方法的优缺点。将小波分析引入干涉纹图噪声抑制当中,研究了基于小波变换和圆周期均值/中值的干涉纹图滤波方法,提出了一种基于小波变换和加权圆周期中值的干涉纹图滤波方法。在综合考虑高频系数方向特性和统计特性的基础上,研究了一种基于小波变换和方向中值的干涉纹图滤波方法,该方法在小波变换的基础上,分别对各个方向的高频系数部分检测其是否为对应方向的边缘,对不同的边缘方向用不同的模板平滑后再进行中值滤波处理。利用德国欧洲空间局(ESA)提供的ERS-1/2串行SAR数据进行了实验,实验结果表明本文所提出的干涉纹图噪声抑制方法不但有效抑制了干涉图中的相位噪声,而且很好地保留了干涉纹图中的边缘细节信息,保证了干涉条纹的平滑连贯,为以后的数据处理和研究工作奠定了基础。关键词:合成孔径雷达干涉;干涉纹图;噪声抑制;小波变换;边缘方向检测II 硕士学位论文AbstractInterferometricsyntheticapertureradar(InSAR)imagingtechniquehasmanyadvantagesoverotherimagingmethods,suchaslargecoverage,highspatialresolution,highelevationaccuracy,andthecapabilityofworkingatalltimeandunderall-weatherconditions.Itisaneffectivemethodtoacquirethethree-dimensionalDigitalElevationModel(DEM).Soitcanbewilelyusedfortopographicsurveying,earthquakemonitoring,volcanomonitoring,deformationdetectionandsoon.BysummarizingthedevelopmenthistoryandthepresentsituationofInSARthree-dimensionalimagingtechnique,analyzingtheprincipleofInSARsystemandmainfactorswhichaffecttheInSARmeasuringaccuracy,thedataprocessingflowofInSARsystemispresented.AccordingtothedifferenceofinterferencephasefromtwoSARpluralimages,theInSARtechniquegainsthethree-dimensionalinformationofgroundobject,thereforethenoiseofinterferogramisoneofthemostimportantfactorswhichaffectthemeasuringaccuracyofthesystem.Inordertoobtainaccuratedataofinterferencephasetomaketwo-dimensionalphaseunwrappingdone,wemustdenoiseeffectively.SotheresearchonInterferogramdenoisingalgorithmofInSARisneeded.Byresearchingtheorigin,characteristicoftheInSARinterferogramnoiseandstatisticalcharacteristicofphase,theadvantageanddisadvantageoftraditionalinterferogramfilteringmethodsareanalyzedbasedonexperimentresults.Byinvolvingthewaveletanalysisintointerferogramdenoising,someresearchworkaboutinterferogramdenoisingmethodofInSARbasedonwavelettransformandperiodicmean/medianalgorithmisaccomplishedandaninterferogramdenoisingmethodbasedonwavelettransformandweightedperiodicpivotingmedianfilteralgorithmisproposed.Inthesynthesisconsiderationonthestatisticalcharacteristicanddirectionalcharacteristicofhighfrequencycoefficient,aninterferogramdenoisingmethodbasedonwavelettransformanddirectionalmedianisproposed.Themethodisbuiltonwavelettransform,thenchecksthehighfrequencycoefficientofalltheorientations,judgeswhetherthecoefficientistheedgeofcorrespondingorientationornot,smoothesthecoefficientusingdifferentmodecorrespondingtoitsorientationandfiltersthecoefficientusingmedianfiltering.III InSAR干涉纹图噪声抑制方法研究ExperimentshavebeencarriedoutbasedonERS-1/2SARtandemdataandtheresultsshowthattheproposedmethodscanreduceinterferogramnoiseefficientlyaswellaspreserveedgeanddetailverywell.Themethodsalsoguaranttheinterferencefringesmoothandestablishthefoundationforthelaterdataprocessingandresearchwork.Keywords:InSAR;Interferogram;Denoising;Wavelettransform;EdgeorientationexaminationIV 湖南大学学位论文原创性声明本人郑重声明:所呈交的论文是本人在导师的指导下独立进行研究所取得的研究成果。除了文中特别加以标注引用的内容外,本论文不包含任何其他个人或集体已经发表或撰写的成果作品。对本文的研究做出重要贡献的个人和集体,均已在文中以明确方式标明。本人完全意识到本声明的法律后果由本人承担。作者签名:日期:年月日学位论文版权使用授权书本学位论文作者完全了解学校有关保留、使用学位论文的规定,同意学校保留并向国家有关部门或机构送交论文的复印件和电子版,允许论文被查阅和借阅。本人授权湖南大学可以将本学位论文的全部或部分内容编入有关数据库进行检索,可以采用影印、缩印或扫描等复制手段保存和汇编本学位论文。本学位论文属于1、保密□,在______年解密后适用本授权书。2、不保密R。(请在以上相应方框内打“√”)作者签名:日期:年月日导师签名:日期:年月日I 硕士学位论文第1章绪论1.1InSAR技术概述[1-3]合成孔径雷达(SAR:SyntheticApertureRadar)是一种主动式高分辨率空间微波遥感成像雷达,它具有全天时、全天候、高分辨率成像的能力,而且作用距离远,成像范围大,能够穿透植被和地表层,解决了光学图像总受天气干扰的难题,因此受到世界各国的广泛重视,已经被成功应用于地质、水文、海洋、测绘、环境监测、农业、林业、气象、军事等领域。随着SAR应用范围的不断扩大,人们除了希望得到目标的几何位置和几何形状等二维信息外,还希望能获取目标的高度信息。但是,传统的SAR技术只能获得目标的二维信息,它缺乏获取地面目标三维信息和监测目标微小形变的能[4-16]力。通过将干涉测量技术与传统SAR技术结合而形成的干涉合成孔径雷达(InSAR:InterferometrySyntheticApertureRadar)三维成像技术提供了获取地面目标三维信息的全新方法,它通过两副天线同时观测(单轨道双天线模式)或通过一副天线两次平行观测(单天线重复轨道模式),获取地面同一景观的复图像对。由于目标与两天线位置的几何关系,在复图像上产生了相位差,形成干涉纹图,干涉纹图中的相位值包含了斜距向上的点与两天线位置之差的精确信息。因此,利用传感器高度、雷达波长、波束视向及天线基线距之间的几何关系,就可以精确的测量出图像上每一点的三维位置和变化信息。InSAR技术可以全天侯大面积地获取高精度的三维地形数据,精度在1m量级,地壳表面变化的测量精度可达1cm。这样高精度的数字地图,对我国的国民经济建设和国防事业具有十分重要的意义,可以广泛应用于地形测量、地质灾害监测、地表沉降监测、海洋表面监测、森林资源调查、土地利用动态监测、军事目标识别和监测等领域。InSAR技术用于获取地面目标三维信息具有许多优势:首先,InSAR技术继承了SAR系统的全天时、全天候和大范围的特点,可以穿透云层或其它遮挡层,不受气候、昼夜等因素影响进行大范围的三维成像,这是一般方法所无法相比的;其次,InSAR技术包含了SAR系统的波长特性,由于SAR系统采用的波长一般为厘米级,因此可以获得很高的测量精度。目前,对于机载InSAR系统可以获得1米量级的测量精度,对于星载InSAR系统精度在5米左右,利用差分合成孔径雷达干涉(D-InSAR:DifferentialInSAR)成像技术对地表形变监测的检测精度可以达到厘米量级;第三,InSAR技术是在SAR系统的基础上发展起来的,1 InSAR干涉纹图噪声抑制方法研究因此它的空间分辨率能达到SAR系统本身所具有的空间高分辨率;第四,InSAR技术可以形成数据处理过程的自动化,从而避免一般光学方法中复杂的人工操作程序。InSAR技术是一个交叉性很强的新领域,涉及对地观测、电磁波传播、信号处理、影像处理与空间大地测量和数字摄影测量等多个领域。因此,数据的处理是一个复杂的过程,从初始的InSAR影像对到提取地面高程信息或地表形变信息的数据处理过程,包含很多环节,而其中的许多关键问题有待进一步的研究与完善。美国航天局、欧洲航天局从90年起连续发射了多颗携带合成孔径雷达的卫星,投入大量人力物力进行InSAR技术的研究。1.2InSAR技术的研究发展20世纪70年代人们就已经提出InSAR的概念,但限于当时的雷达和计算机技术水平,无法大量获得相干性好的干涉数据并进行有效的处理。近年来,随着上述领域科技的突飞猛进,InSAR己经从理论研究阶段进入应用阶段,并在地形测绘、地震和火山监测以及地形沉降检测等方面取得大量富有开拓性意义的成果。[13]干涉雷达最早是在1969年应用于金星和月球表面状况的探测。1974年,[6]Graham第一次提出应用合成孔径雷达进行地形制图,地表的三维数据可以通过机载或星载的SAR传感器干涉测量的几何特征测试来获得,但下一篇公开发表的论文直到1985年才由美国加利福尼亚州的喷气推进实验室(JPL:JetPropulsion[11]Laboratory)提出。到1986年,Zebker和Goldstein利用机载侧视雷达获得了第一个使用的观测结果。他们在飞机上安装了两根相距11.1米的SAR天线(两根天线间的距离称为基线(Baseline))进行干涉测量,由一部天线发射微波信号,第二部天线同步接收回波信号,在未校正飞机状态畸变的情况下,得到的复图像的分辨率大约为10米左右。测得的地表高程数据在海洋上其均方误差在2-10米左右,而在陆地上的一些区域其误差则更大一些。Goldstein等人在1988年将星载雷达的干涉测量思想应用于SEASAT卫星,对美国死谷格登堡盆地地区进行了实验研究,得到的地表地形数据与美国地质调查局出版的勘测图吻合得相当好。1990年Li和Goldstein利用星载SAR的重复轨道数据进行了干涉测量实验,并研究了基线距对干涉结果的影响。1991年欧洲空间局(ESA:EuropeanSpaceAgency)发射了ERS-1卫星,该卫星搭载了一个C波段的合成孔径雷达,这成为InSAR研究的一个重要的里程碑。由于有了优质易得的InSAR数据源,大批欧洲研究者加入到这个领域,亚洲(主要是日本)的一些研究者也开展了这方面的研究。大量的文章探讨了InSAR应用潜力和制约因素,并进行了各种模拟实验来改善由来的研究。到1995年,ERS-2 硕士学位论文2卫星发射上天,与ERS-1配合,特别安排了双星追逐观测期(tandemmission),为InSAR的研究提供了数据保证。利用ERS-1和ERS-2两颗卫星的重复轨道观测,可以获得1天间隔的SAR数据,这极大的推动了星载SAR的干涉测量研究。一直到现在,这些测量数据的研究和结果确认工作仍在进行,并且取得了巨大的进展。美国奋进号航天飞机分别于1994年4月和10月两次升空,执行由美国、德国、意大利等13个国家参加的航天飞机成像雷达计划(SIR-C/SAR)。该计划成功地获取了大量适合干涉测量的多波段全极化SAR数据,进行了干涉测量的多波段全极化研究,以及在地形、火山检测中的应用研究。在上世纪90年代的InSAR研究热潮中,可以这样总结:在星载SAR研究方面,SEASAT、ERS-1/2、JERS-1、SIR-CIX-SAR、RADARSAT等数据都进行过重复轨道干涉测量实验;在机载系统研制方面,目前已有的机载干涉雷达系统有美国NASA/JPL的AIRSAR/TOPSAR系统、丹麦的EMI-SAR系统,它们都装载两部天线,一次飞行即可获得高程信息,测高精度0.1-1米。目前可用于InSAR研究的数据以ERS-1/2为佳,其次是SIR-C/X-SAR和RADARSAT,JPL的TOPSAR是机载系统的佼佼者。2000年2月11日至22日美国所执行的SRTM(ShuttleRadarTopographyMission)全球地貌测绘计划是InSAR研究应用领域中又一个重要的里程碑。该计划通过航天飞机环绕地球11天的飞行,利用C波段及X波段SAR获得了覆盖地球表面80%陆地范围的地貌数据,经进一步处理后,至少可以得到绝对精度在10-15m左右的、几进全球范围的数字高程地图。此次测绘的结果不但具有重要的科学及经济价值,对未来军事活动的作用也不可低估。[28-32]国内从1996年开始InSAR方面的研究,主要开展了一些卫星雷达干涉测量原理、成像技术、相位解缠、基线参数估计及地形高度估计等研究,提出了一些新的方法,并在一些应用领域(地形测量、城市变化监测、获取数字地形高度DEM图)取得成功。近年来,我国通过参加SIR-C/X-SAR计划、ERS计划、RADARSAT、ADRO计划等国际星载SAR计划,开展了干涉雷达数据处理和应用研究,利用航天飞机成像雷达数据,提取了试验区数字模型,利用中科院电子所研制的机载SAR系统,通过了干涉测量试验,并且已有一些论文在国内外重要刊物上发表。但从目前研究水平来看,我国距世界发达国家研究水平仍有很大差距,还处在刚起步的阶段,还远远没有形成系统,更谈不上广泛应用,需要尽快研制我国星载SAR,同时利用国外星、机载InSAR数据,开展InSAR处理算法和应用的研究,使我国干涉雷达的研究和发展能够在国际上占据应有的位置,并更好地服务于国民经济和国防建设。3 InSAR干涉纹图噪声抑制方法研究1.3InSAR技术的应用InSAR成像技术利用SAR复图像中含有的相位信息,通过干涉处理来提取目标的三维信息,因此用于制作地形图、生成高精度的数字高程图DEM是自InSAR技术研究和应用以来的主要应用领域。D-InSAR成像技术由于能获得厘米甚至毫米级的高精度三维形变,因此可应用于DEM的修测与精化、地震变形、火山变形、地面沉降、滑坡、冰川移动、军事目标识别等方面。1.3.1DEM的建立和优化对地表面形态的研究,即地形、地貌的研究,是一个最基础、最重要的研究,而获得高质量的DEM是其关键。InSAR技术用于建立DEM的优势在于它能提供比其它手段更高的精度和更快的速度。InSAR技术的测量精度由于成像几何和干涉图像质量不同而有较大的波动,精度好的可达米级。研究结果表明,InSAR技术用于获取DEM是非常有效的,特别是在人烟稀少、环境恶劣和光学遥感图像难以获得的地区(如热带雨林、火山地区、极地地区等),InSAR技术更是一种有效的测绘手段。利用InSAR快速准确地获得全世界的数字高程图,这在世界范围内都是具有重要意义的,构建数字化地球的梦想将不再遥远。同时,利用InSAR技术获取的高精度DEM可以对原有的低精度DEM进行优化,通过与其它技术获取的DEM进行融合,既可以使整体的DEM精度得到提高,还可以弥补由于雷达阴影和雷达折叠的影响造成InSAR技术无法生成某些区域DEM的缺陷。1.3.2地质灾害研究地震、火山喷发等突发性的地质灾害,严重影响着人类的生命和财产安全,研究这些地质灾害的成因及规律已成为地球科学研究者的重要课题。地震灾害是地壳形变积累的结果,对地壳形变过程的监测和震后测量具有重要的意义。目前,监测地震灾害主要依赖于VLBI技术、SLR技术、GPS技术、水准测量、应力测量和张力测量等传统方法,但这些方法局限于有选择性地监测一些离散站点。D-InSAR技术则具有高分辨率、高精度、连续空间覆盖等特点,它通过对卫星获得的震前和震后的SAR复图像进行干涉处理,可以提供有关地震的震源位置、震源机制和地震滑动分布等信息,为地震研究提供可靠和宝贵的数据。InSAR技术也可用于检测与火山喷发有关的地面形变模式。从由InSAR技术获得的火山的高分辨率DEM,可以估计火山熔岩流场的结构和形态的细节。它不仅可以作为这些火山将来由于喷发而引起变化的检测基准,而且也可用于估计未来的熔岩流、泥流和火山碎屑流所引起的灾害。4 硕士学位论文地质构造、板块运动和人为因素(如城市发展、地下水抽取等)都会导致地面的沉降,城市地面沉降问题一直受到各国政府的高度重视。采用D-InSAR测量技术可以在大面积范围内(100km×100km)监测地面的微小形变,具有不需要人员进入现场区域测量的特点,相对于GPS测量和水准测量手段而言,有很大的成本优势,具有广阔的市场应用潜力。1.3.3用于极地监测极地冰盖和冰川对地球气候的变化起着极其重要的作用,因此对极地冰盖体积和冰川运动及其边缘位置变化的监测具有非常重要的意义。冰川的运动速度相对较高(一般平均每天几十厘米),因此,监测时间间隔不能太长。与传统的监测手段相比,D-InSAR技术具有大范围、高效率等特点,它可以精确快速地测量极地冰盖厚度的变化和冰川的移动情况。1.3.4用于海洋监测海洋占据地球表面的70%,蕴藏着人类赖以生存的重要资源。而海面上的天气状况往往非常恶劣,给光学遥感手段监测海况带来极大困难。SAR是进行海洋观察的最理想工具之一,通过InSAR技术,不仅可以探测到海面上船舶的运动方向和速度,而且可以观察到多种海洋动力学现象以及由于海底地形变化而引起的海面波浪的差异等。总之,经过十几年来的发展,InSAR成像技术和D-InSAR成像技术的研究与应用正在不断地完善与深入。我们相信,随着这项新技术的不断发展,它必将给人类带来巨大的经济和社会效益。1.4InSAR的未来随着ERS-1、ERS-2以及SIR-C等卫星的成功发射并获取了比较理想的数据,以航空航天InSAR为主的InSAR技术蓬勃的发展起来。2000年,SRTM成功的获取了覆盖地球陆地绝大部分的InSAR数据,SRTM被认为是第一颗真正的InSAR卫星,是InSAR发展的一个里程碑。InSAR的未来一片光明,新的TOPSAT卫星计划已在孕育之中。JPL已为提出的TOPSAT计划做了许多准备工作,包括在机载TOPSAT上的预演试验,如GPS-INS导航定位。TOPSAT将是一个双星系统,但不同于ERS-1/2,两颗星不是间隔一天,而是只相隔1km,预计其轨道设计将更加优化。TOPSAT将有多波束激光测高仪作为补充,卫星设计寿命为3年,两颗卫星可在3个月内获取南纬70度到北纬70度的InSAR数据。两极地区由激光测高仪补测,InSAR基线随纬度不同而异,在赤道为2km,随着纬度增大,基线逐渐减少,最小为700m,星上5 InSAR干涉纹图噪声抑制方法研究将载有GPS-INS导航定位和姿态定向系统,以精确确定卫星位置、姿态和InSAR基线。GPS定位精度将在10cm,InSAR数据的精度将达到平面精度30m,高度精度2~5m,激光测高仪在起伏不大的地区预计可达到1m精度,计划选择某些地区由激光测高仪量测植被高度和地表粗糙度。TOPSAT的设计和不久将来的升空运行,是InSAR发展历程中的又一个重大事件,对全球环境、资源、经济社会可持续发展,灾害预报、防治和监测都将具有重要的意义。1.5论文的结构安排本论文的研究工作在国家自然科学基金(60375001)和教育部重大项目培育基金(706043)的资助下,首先根据InSAR三维成像技术的特点,分析了InSAR系统的工作原理和引起InSAR数据测量误差的各种因素,然后从InSAR系统的数据处理流程出发,对InSAR干涉相位纹图的噪声抑制方法进行了重点研究。在研究InSAR干涉纹图噪声来源和噪声特征的基础上,分析了传统的干涉纹图滤波方法的优缺点,提出了基于小波分析的干涉纹图噪声抑制方法,实验结果表明所提方法不仅能有效地抑制干涉纹图中的噪声,而且能很好地保留了干涉纹图中的边缘细节信息,为以后的数据处理和研究工作奠定了基础。论文主要由以下几部分组成:第1章从研究课题的背景入手,介绍了InSAR技术的研究意义及国内外现状。第2章详细介绍InSAR技术的工作模式、基本原理、数据处理流程以及影响InSAR测量精度的各种因素。第3章对InSAR干涉纹图中的噪声进行了分析,包括其来源、特性、相位统计特征等,并介绍了评价干涉纹图质量的标准。第4章总结了空间域的干涉纹图噪声抑制方法,以实验结果为基础,全面分析了这些算法的优缺点。第5章将小波分析方法引入InSAR干涉纹图噪声抑制,研究了基于小波变换和圆周期均值/中值、基于小波变换和加权圆周期中值、基于小波变换和方向中值的干涉纹图滤波方法,实验结果验证了所提方法的有效性和优越性。结论对全文的工作进行了总结,同时对今后进一步的工作提出了一些建议。6 硕士学位论文第2章InSAR技术原理2.1InSAR的工作模式根据InSAR平台和使用条件的不同,InSAR的工作模式可以分为交叉轨道干涉(across-trackinterferometry)测量模式、沿轨道干涉(along-trackinterferometry)测量模式和重复轨道干涉(repeat-trackinterfermetry)测量模式三种。2.1.1交叉轨道干涉测量模式交叉轨道干涉测量模式是飞行平台上同时装载两个天线,其中一个负责发[13]射并接收雷达波束,另一个则负责接收,这样一次飞行即可获得干涉影像对。在基线固定,能准确确定平台位置的情况下,就能获得高质量的干涉测量数据和高程计算结果,在航空平台上多采用这种方式。它的优势在于精度高而且机动性能好。其干涉几何原理图如图2.1所示,从图中可以看出,两副天线的安装位置与飞行方向垂直。在该模式下,干涉相位差是由于地面目标的高度变化引起的,所以主要用于地形制图和地形变化监测。z飞行轨道Hθrr21xhy图2.1交叉轨道干涉测量模式2.1.2沿轨道干涉测量模式沿轨道干涉模式与交叉轨道干涉模式一样,都要求两副天线安装在同一平台上,因此目前也主要适用于机载SAR系统。其干涉几何原理图如图2.2所示,此时两副天线沿飞行方向相隔一段距离。采用该模式得到的相应像素的相位差是由于测量时物体的运动产生的,因此它适用于对运动的目标进行监测,如海洋制图、波浪谱测量等。7 InSAR干涉纹图噪声抑制方法研究zHθr飞行轨道r21xhy图2.2沿轨道干涉测量模式2.1.3重复轨道干涉测量模式重复轨道干涉测量模式只要求一根天线,在尽可能短的时间间隔内,经过几乎相同的轨道,以微小的几何视角差两次获取同一地区的数据,由于受大气的影响较小,卫星比飞机具有更准确、更稳定的飞行轨道,因此这种方式最适合800千米以上高空的卫星。目前ERS-1/2双星重轨模式因只相隔一天,成为采用这种模式最成功的卫星。该模式的优势在于能够快速获取大范围或全球范围的干涉数据。目前采用此模式的还有日本JERS-1装载的SAR系统和美国航天飞机上的SIR-C/X-SAR系统,也都取得了很好的效果,其干涉几何原理图如图2.3所示。oz2飞行轨道o1rHθ2r1xhy图2.3重复轨道干涉测量模式2.2InSAR三维成像技术原理InSAR三维成像技术是一门根据复雷达图像的相位数据来提取地面目标三维空间信息的技术,其基本思想是:利用两副天线同时成像或一副天线相隔一定时8 硕士学位论文间重复成像,获取同一区域的复雷达图像对,由于两副天线与地面某一目标之间的距离不等,使得在复雷达图像对同名像点之间产生相位差,形成干涉图像,干涉图像中的相位值即为两次成像的相位差测量值,根据两次成像相位差与地面目标的三维空间位置之间存在的几何关系,利用飞行轨道的参数,即可测定地面目标的三维坐标,它可以用来提供大范围的高精度数字高程模型(DEM)。本文以卫星重复轨道干涉模式为例,其成像几何示意图如图2.4所示。S2BSα1BB⊥||r+δrθrHPh图2.4INSAR的成像几何示意图假设S1和S2是卫星两次对同一地区成像的天线位置,两次位置之间的距离B称为基线距,P为地面上高度为h的一目标,V1(r,x)和V2(r,x)分别为点目标P在S1和S2成像后的复图像中的复值,则可得:V(r,x)=V(r,x)ejφ1(2.1)11V(r,x)=V(r,x)ejφ2(2.2)22式中φ1和φ2分别为P点在S1和S2的回波相位。将在S1和S2成像的两幅复图像进行干涉得:V(r,x)⋅V*(r,x)=V(r,x)V(r,x)ej(φ1−φ2)=V(r,x)V(r,x)ejΦ(2.3)121212式中*表示复共轭,Φ为P点在S1和S2中的回波相位差,即干涉相位。由式(2.3)可得:*Φ=arctan[V1(r,x)⋅V2(r,x)](2.4)设S1和S2到P点的路径长度分别为r和r+δr。则P点在两幅SAR复图像中的相位分别为:4πrφ1=(2.5)λ9 InSAR干涉纹图噪声抑制方法研究4π(r+δr)φ2=(2.6)λ式中λ为雷达信号的波长。由式(2.5)和(2.6)可得两次测量的相位差Φ为:4πΦ=φ1−φ2=−δr(2.7)λ上式体现了两次雷达成像的相位差与雷达信号到地面目标信号路径差的关系,而等式左边的两次测量相位差Φ可以通过式(2.3)由两幅复图像生成的干涉图像来获得。假设视角为θ,基线距与水平方向的夹角为α,根据余弦定理可得:222π22(r+δr)=r+B−2rBcos(−θ+α)=r+B−2rBsin(θ−α)(2.8)22因δr<0(5.1)x∫a−∞a等效的频域表示是:a+∞*jωτWTa(,)τ=X()()ωψaωedω(5.2)x∫2π−∞式中,X()ω,ψ()ω分别是x(t),()ψt的傅立叶变换。小波变换的时频窗口特性与短时傅立叶的时频窗口不一样,因为τ仅仅影响窗口在相平面时间轴上的位置,也影响窗口的形状。这样小波变换对不同的频率在时域上的取样步长是可调节的,即在低频时小波变换的时间分配率较低,而频率分辨率较高;在高频时小波变换的时间分辨率较高,而频率分辨率较低,这正符合低频信号变化缓慢而高频信号变化迅速的特点。这便是它优于经典的傅立叶变换和短时傅立叶变换的地方,从总体上说,小波变换比短时傅立叶具有更好的35 InSAR干涉纹图噪声抑制方法研究时频窗口特性。[42]由此可见,小波变换具有以下特点和作用:(1)具有多分辨率(也叫多尺度)的特点,可以由粗到细地逐步观察信号。(2)我们也可以把小波变换看成用基本频率特性为ψ()ω的带通滤波器在不同尺度a下对信号做滤波。由于傅立叶变换的尺度特性,如果ψ()t的傅立叶变换是tψ()ω,则ψ()的傅立叶变换为aaψ()ω,因此这组滤波器具有品质因数恒定,a即相对带宽(带宽与中心频率之比)恒定的特点。(3)适当地选择小波基,使ψ()t在时域上为有限支撑,ψ()ω在频域上也比较集中,便可以使小波变换在时、频两域都具有表征信号局部特性的能力,这样就有利于检测信号的瞬态或奇异点。5.1.2小波应用于图像的分解与重构图像是二维信号,因而将小波引入图像的分解和一维信号的分解是不同的。一幅图像,其中任意一点(x,y)有一个图像信号的灰度值f(x,y)与之对应,点坐标(x,y)连续变化时就确定了一个连续变化的二维信号函数f(x,y)。要将小波变换应用到图像处理中,就必须首先把小波变换由一维推广到二维,即对f(x,y)进行信号处理。二维小波的图像分解我们可以用如下的论述来推导。11我们假设一个一元尺度函数φ(x)生成一个多分辨分析{V}而一元尺度函数k2212φ(y)生成另一个多分辨分析{V},则V与V的张量积空间kkk12VVV=⊗(5.3)kkk1kk/212kk/22由于V的底基是{2(φ2)x−j},而V的底基是{2(φ2)yl−},所以V的底基kkkkk/212k是{2φφ(2x−−jy)(2l)}。111222设V关于V的补空间是W,V关于V的补空间是W,即kk+1kkk+1kg111VVW=+(5.4)kkk+1g222VVW=+(5.5)kkk+1通过计算有:gVVW=+(5.6)kkk+1gg(1)(2)(3)WVW(1)=⊗12WWV(2)=⊗12其中,WWWWkkkk=++,而kkk,kkk,36 硕士学位论文(3)12WWWkkk=⊗2任给fxyL(,)∈(IR),(,)fxy是f在V中的投影。这时,对f(,)xyV∈,NNkkgxyW(,)∈,有kkf(,)xy=fxygxy(,)+(,)(5.7)kk+1k而gxy(,)还可进一步分解为k(1)(2)(3)gggg=++(5.8)kkkk()ii()其中gW∈(i=1,2,3)。kki设{}a{}b(i=1,2,3)是由两个一元分解序列生成的二元分解序列lj,lj,12112aa==a,bablj,llj,jlj(5.9)212312ba==a,bbblj,llj,jlj且kkfxykk(,)=−∑c;,nmφ(2xnym,2−)nm,(5.10)iiikkgxykk(,)=−∑d;,nmψ(2xnymi,2−),=1,2,3nm,则有分解算法:cacknm;,=∑lnjmklj−−+2,21;,lj,(5.11)iidbcknm;,=∑lnjmklj−−+2,21;,lj,图5.1是图像分解的示意图,L表示低频,H表示高频,下标1,2表示一级或二级分解。[44]小波的分解以及抽样的示意图如图5.2所示。图5.2中,2↓1表示对列下抽样,保留列号为偶数的值;1↓2表示对行下抽样,保留行号为偶数的值。从图5.1我们看到,原图像经过一级小波分解,分解为LL1、LH1、HL1、HH1四个分量,其中LH1、HL1、HH1分别为水平方向、垂直方向、斜线方向的高频分量,LL1为可以进行再分解的低频分量,小波的二级分解就是将LL1分解为LL2、LH2、HL2、HH2四个分量,如果再进行第三级分解则是将LL2按照同样的方法再进行分解。i同样,设{}p,{}q(i=1,2,3)是由两个一元两尺度序列得到的二元两尺度序lj,lj,列,即:12112p==pp,qpqlj,llj,jlj(5.12)212312q==qp,qqqlj,llj,jlj37 InSAR干涉纹图噪声抑制方法研究LL1HL1原图LH象LH1HH1LL2HL2HL1LH2HH2LH1HH1图5.1图像小波分解示意图行卷积列卷积C(n,m)LL121am()21an()12LH11bn()12HL11bm2()21an()12HH1112bn()图5.2小波分解及其下抽样[44]则有重构算法:3iicpkn+−1;,m=+∑∑()nl2,mjk−−2c;,ljqnl2,mjk−2d;,lj(5.13)lj,1i=图5.3是小波抽样及重构的示意图:38 硕士学位论文行卷积列卷积LL1C(n,m)112pm2()21pn()LH1121qn()HL121112qm2()pn()HH1121qn()图5.3小波上抽样及其重构图中,2↑1表示对列上抽样,奇数列号位置插入0;1↑2表示对行上抽样,奇数行号位置插入0;图5.2把尺度为j的c(n,m)的低频部分分解成尺度为j+1的低频部分LL1和水平方向、垂直方向、斜线方向的3个高频部分,再依此将LL1分解。同理,小波多分辨分析的重构将尺度为j+1的低频部分和3个方向上的高频部分的系数进行处理,从而得到重构后的图像。实验采用德国欧洲空间局(ESA)提供的我国新疆乌鲁木齐市附近地区的ERS-1/2串行SAR数据,两幅单视SAR复图像经过配准后进行干涉处理得到干涉纹图,取大小为512×512像素的子图像作为实验图进行小波的分解与重构,实验结果如图5.4所示(以二级分解为例):a)原始干涉纹图39 InSAR干涉纹图噪声抑制方法研究b)第一级分解后得到的四个方向上的分量图c)第一级分解后得到的低频分量LL140 硕士学位论文d)对LL1进行第二级分解后得到的四个方向上的分量图e)第二级分解后得到的低频分量LL2图5.4乌鲁木齐地区干涉纹图的二级小波分解从图5.4不难看出,对干涉相位图进行到第二级小波分解时,高频分量里的噪声斑点仍然清晰可见,而其低频分量里的噪声斑点明显少于高频分量里的噪声斑点。因此,干涉相位图的噪声主要存在于高频分量之中。图5.5为二级分解后的干涉相位图的二级重构以及与原图的对比。41 InSAR干涉纹图噪声抑制方法研究a)一级分解后的LL1原图b)二级重构出的LL1图c)原始干涉纹图d)一级重构出的原始干涉相位图图5.5干涉纹图的二级小波重构从图中可以看出,重构出的图像与原图的差别很小,小波变换是一种有效的分解图像的高低频并重构的方法。5.1.3噪声在小波分解下的特性设2维图像信号模型可以表述为:s(i,j)=f(i,j)+σge(i,j),i,j=0,,m-1ggg式中,s(i,j)为真实信号,e(i,j)为噪声信号,是标准偏差不变的高斯白噪声。s(i,j)为含有噪声的图像信号。将噪声看成是一个普通的信号,假设对噪声进行小波分解的系数用c(j,k)表示,j表示小波尺度,k表示时间。因此,噪声信号具有以下特性:1.如果所分析的信号是一个平稳、零均值的白噪声,则其小波分解系数是不相关的。2.如果所分析的信号是一个高斯噪声,则其小波分解系数是独立的,且也是高斯分布。3.如果所分析的信号是一个有色、平稳、零均值的高斯噪声序列,则其小波分解系数也是高斯序列。对每一个分解尺度j,其系数是一个有色、平稳的序列。4.如果所分析的信号是一个固定的零均值ARMA模型,则对每一个小波分解42 硕士学位论文尺度j,其小波分解系数也是固定的、零均值的ARMA模型,其特性取决于j。5.2基于小波变换和圆周期均值的滤波方法通过前面对小波的介绍,可以看出小波分析能同时在时域和频域中对信号进行分析:在频域内分辨率高时,时域内分辨率则低;在频域内分辨率低时,时域内分辨率则高,具有自动变焦功能。本文所研究的对象是消除了平地效应后的相位图,其中存在大量的噪声信号。设干涉相位图的实部和虚部两个数据集中的噪声为标准差不变的高斯白噪声。针对噪声在小波分析下的特性,由于小波分析具有能有效的区分信号中的突变部分和噪声的特点,因此,采用小波分析的方法对干涉相位图进行滤波既能够滤去噪声,又能够保持相位突变时的相位值。应用乘性噪声模型X=SNg(5.14)式中,S表示信号,为一高斯随机过程,均值和方差分别为µs和σs;N表示噪声,是一均值为1的高斯随机过程,且两者是相互独立的随机过程。X表示原始数据。由公式(5.14)有µµ=XS22(5.15)2KKXN−K=S21+KN式中,KKK===σ///µσµσµ为归一化标准差。因此,可以利SSSNNNSSS,,用K的局部估计值来区分K是否为均匀区域或非均匀区域。XS基于小波分析的这些特性,本文提出了基于小波变换和圆周期均值算法的噪声抑制方法。5.2.1窗口的选择选择窗口的大小实际上与系统的参数有关。根据前文所述,我们知道干涉纹图存在的几个主要噪声:热噪声为InSAR系统本身固有的噪声来源,通过改进SAR系统本身的设计可以降低噪声影响;相干斑点噪声随机性很大,大大降低图像的信噪比,它的大小受目标表面相对于系统波长的粗糙程度所影响;基线失相关(Spatialdecorrelation)和时间失相关(Timedecorrelation)对噪声大小的影响在于:ò基线长度为零时,则不存在基线失相关,则在此理想状态下,有基线引起的噪声程度为0;当基线长度到达某一临界值时,两次雷达回波完全没有相关性,导致噪声程度非常高。ò重复轨道的两次成像时间不同,会导致地面目标的散射特性随环境的变化而43 InSAR干涉纹图噪声抑制方法研究变化,时间越长,则时间失相关程度越强,这意味着,时间越长,两幅数据图像产生出的干涉纹图所含噪声越高;反之亦然。一般来说,重复轨道的观测周期越短,由时间失相关引起的噪声就越小,为了减小时间失相关噪声的影响,应尽量选取时间间隔小的两幅复图像数据。本文的实验数据是ERS-1/2串行采集得到的一对复数图象数据,两幅图像数据采集的间隔时间为1天,即时间失相关很小,但基线距离相对较大。为了比较不同大小窗口对滤波效果的影响,本文采用了固定窗口(窗口大小为7×7)和i−1变化窗口(窗口大小为ww=−21(2i≥),其中w=7为第一层小波分解后对高i00频分量进行均值滤波的窗口大小,i为小波分解的层数)两种窗口进行实验对比。5.2.2算法的实现步骤基于小波变换和圆周期均值的滤波方法的实现原理如图5.6所示:低频分量LLLL1LL2原始干涉小波小波相位图滤波后的干分解涉相位图重构HL,LH,LL圆周期均值滤波后的HL1,LH1,HH1高频分量HL2,LH2,HH2高频分量图5.6基于小波变换和圆周期均值的滤波方法原理图算法具体的实现步骤如下:1.先对干涉纹图进行小波分解。选择一个小波(本文选择sym4小波)和小波分解的层次N(本文N=3),然后计算分解后各层的系数。2.对各层的各高频分量分别进行圆周期均值滤波,采用第4章公式(4.4)的方44 硕士学位论文法实现。窗口的选取采用上述的固定窗口和变化窗口两种方式。3.干涉纹图的小波重构。根据小波分解后的第N层的低频系数和经过滤波后的第N层的各高频系数重构出第N-1层的低频系数,再依此类推来进行2维信号的小波重构,重构所得图像即为滤波后的干涉相位图。5.2.3实验结果与分析基于小波变换和圆周期均值的滤波方法实验结果如图5.7和5.8所示。a)喀什地区原始干涉纹图b)7×7固定窗口小波和圆周期均值滤波后的图像c)变化窗口小波和圆周期均值滤波后的图像图5.7喀什地区小波和圆周期均值滤波效果图a)乌鲁木齐地区原始干涉纹图45 InSAR干涉纹图噪声抑制方法研究b)7×7固定窗口小波和圆周期均值滤波后的图像c)变化窗口小波和圆周期均值滤波后的图像图5.8乌鲁木齐地区小波和圆周期均值滤波效果图实验前后干涉纹图中的残差点比较如表5.1所示。表5.1小波和圆周期均值滤波残差点比较表噪声抑制方法残差点数喀什地区原始干涉纹图(无滤波)54428小波和圆周期均值(7×7)1472小波和圆周期均值(一层7×7,2层13×13,三层27×27)1350乌鲁木齐地区原始干涉纹图(无滤波)79684小波和圆周期均值(7×7)2523小波和圆周期均值(一层7×7,二层13×13,三层27×27)2402为了更直观的说明滤波效果,图5.9和5.10所示为滤波前后干涉纹图中残差点的分布比较图。a)喀什地区原始干涉图残差点分布图46 硕士学位论文b)7×7小波和圆周期均值滤波后的残差点分布图图5.9喀什地区残差点分布图的对比a)乌鲁木齐地区原始干涉纹图的残差点分布图47 InSAR干涉纹图噪声抑制方法研究b)7×7小波和圆周期均值滤波后的残差点分布图图5.10乌鲁木齐地区残差点分布图的对比由实验结果可知,基于小波变换和圆周期均值的滤波方法对噪声的抑制能力要优于圆周期均值滤波方法,并且采用此方法滤波后的干涉纹图在相位突变处的平滑明显减少,很好的保留了干涉纹图中的边缘细节信息。但在实验过程中,该方法的计算工作量要大于圆周期滤波方法。5.3基于小波变换和圆周期中值的滤波方法5.3.1算法的实现步骤基于小波变换和圆周期中值的滤波方法,具体的实现步骤为:1.先对干涉纹图进行小波分解。选择一个小波(本文选择sym4小波)和小波分解的层次N(本文N=3),然后计算分解后各层的系数。2.对各层的各高频分量分别进行圆周期中值滤波,采用第4章公式(4.5)的方法实现。本文采用固定窗口(窗口大小为7×7)和变化窗口(窗口大小为i−1ww=−21(2i≥),其中w=7为第一层小波分解后对高频分量进行均值滤波的i00窗口大小,i为小波分解的层数)两种窗口进行实验对比。3.干涉纹图的小波重构。根据小波分解后的第N层的低频系数和经过滤波后的第N层的各高频系数重构出第N-1层的低频系数,再依此类推来进行2维信号的小波重构,重构所得图像即为滤波后的干涉相位图。5.3.2实验结果与分析基于小波变换和圆周期中值的滤波方法的实验结果如图5.11和5.12所示。48 硕士学位论文a)喀什地区原始干涉纹图i−1b)7×7窗口小波和圆周期中值滤波后的图像c)ww=21−≥(2i)窗口小波和i0圆周期中值滤波后的图像图5.11喀什地区小波和圆周期中值滤波效果图a)乌鲁木齐地区原始干涉纹图49 InSAR干涉纹图噪声抑制方法研究i−1b)7×7窗口小波和圆周期中值滤波后的图像c)ww=21−≥(2i)窗口小波和i0圆周期中值滤波后的图像图5.12乌鲁木齐地区小波和圆周期中值滤波效果图实验前后干涉纹图中的残差点比较如表5.2所示。表5.2小波和圆周期中值滤波残差点比较表噪声抑制方法残差点数喀什地区原始干涉纹图(无滤波)54428小波和圆周期中值(7×7)1535小波和圆周期中值(一层7×7,2层13×13,三层27×27)1424乌鲁木齐地区原始干涉纹图(无滤波)79684小波和圆周期中值(7×7)2603小波和圆周期中值(一层7×7,2层13×13,三层27×27)2461由实验结果可看出,基于小波变换和圆周期中值的滤波方法很好的滤除了干涉纹图中的噪声,同时继承了圆周期中值的滤波特点,保留了其在保持相位条纹连续性方面较强的能力。但随着滤波窗口的增长,此方法的计算量急剧增大,计算时间因而变得很长。5.4基于小波变换和加权圆周期中值的滤波方法5.4.1加权圆周期中值滤波算法[48-49]加权圆周期中值滤波算法是在结合已有的一些滤波方法的基础上,利用圆周期均值滤波和圆周期中值滤波算法在滤除除高斯噪声和颗粒噪声等噪声方面的各有所长,将两种算法结合起来,以期望能够同时滤除多种噪声,从而达到较好的滤波效果,其具体的思路如下:(1)在每一个滤波窗口内求出各相位矢量的和矢量d,然后取出d相位值lm,lm,的主值。'(2)对于滤波窗口区域内所有点相角差arg[exp(jQ)/d],首先求其中值,klkm,lm,然后对区域内每一点相角差都以这个中值为基础计算权值,越是接近该中值的相角差权值越大,而且权值归一化,然后将区域内每点相角差与此权值相乘再求和。50 硕士学位论文上述加权圆周期中值滤波的算法实现过程为:lL++DDmM−Qlmout(,)=+∑∑WQklkm,,.klkm,arg(dlm)(5.16)kllLkmmM=−DD=−−1/(1+−QM)klkm,lm,其中W=klkm,Slm,_'Qj=arg[exp(Q)/d]klkm,klkm,,lm'M=median[arg(jQ)/d]lm,,klkmlm,lL++DDmM_2SQlm,,=+∑∑[1/(1(kLkm,−Mlm))]kllLkmmM=−DD=−其余的字母和符号同公式(4.4)。该思路的算法优点在于,当以中值为基础计算权值时,给颗粒噪声点赋的权值非常小,以至于在做累加的时候颗粒噪声的值可以被忽略,这样可以消除一部分颗粒噪声.同时,做累加类似于均值滤波,可以抑制一部分高斯噪声。5.4.2基于小波变换和加权圆周期中值算法的实现步骤基于小波变换和加权圆周期中值的滤波方法的具体实现步骤如下:1.先对干涉纹图进行小波分解。选择一个小波(本文选择sym4小波)和小波分解的层次N(本文N=3),然后计算分解后各层的系数。2.对各层的各高频分量分别进行加权圆周期中值滤波,采用公式(5.16)的方法实现。本文采用固定窗口(窗口大小为7×7)和变化窗口(窗口大小为i−1ww=−21(2i≥),其中w=7为第一层小波分解后对高频分量进行均值滤波的i00窗口大小,i为小波分解的层数)两种窗口进行实验对比。3.干涉纹图的小波重构。根据小波分解后的第N层的低频系数和经过滤波后的第N层的各高频系数重构出第N-1层的低频系数,再依此类推来进行2维信号的小波重构,重构所得图像即为滤波后的干涉相位图。5.4.3实验结果与分析基于小波变换和加权圆周期中值滤波方法的实验结果如图5.13和5.14所示。51 InSAR干涉纹图噪声抑制方法研究a)喀什地区原始干涉纹图i−1b)7×7窗口小波和加权圆周期中值滤波后的图像c)ww=21−≥(2i)窗口小波和i0加权圆周期中值滤波后的图像图5.13喀什地区小波和加权圆周期中值滤波效果图a)乌鲁木齐地区原始干涉纹图52 硕士学位论文i−1b)7×7窗口小波和加权圆周期中值滤波后的图像c)ww=21−≥(2i)窗口小波和i0加权圆周期中值滤波后的图像图5.14乌鲁木齐地区小波和加权圆周期中值滤波效果图残差点比较表如表5.3所示。表5.3小波和加权圆周期中值滤波残差点比较表噪声抑制方法残差点数喀什地区原始干涉纹图(无滤波)54428小波和加权圆周期中值(7×7)1392小波和加权圆周期中值(一层7×7,2层13×13,三层27×27)1286乌鲁木齐地区原始干涉纹图(无滤波)79684小波和加权圆周期中值(7×7)2415小波和加权圆周期中值(一层7×7,2层13×13,三层27×27)2350加权圆周期中值算法集合了圆周期中值和圆周期均值的优点,将它与小波变换相结合的算法在干涉纹图噪声抑制方面取得了很好的效果,残差点的数目进一步减少,同时在保持相位信息和边缘信息上有不错的表现,它的不足在于计算量大,计算时间长。5.5基于小波变换和方向中值的滤波方法5.5.1算法实现步骤基于小波变换和方向中值的滤波方法在考虑高频系数的统计特性的同时,检测各个方向的高频系数的方向特征,通过检测小波变换得到各个方向高频系数的对应方向边缘后,对边缘方向用模板平滑后再中值滤波,因而在抑制噪声的同时保持图像的边缘信息。具体处理步骤如下:1.先对干涉纹图进行小波分解。选择一个小波和小波分解的层次N,然后计算hvd分解后各层的系数,得到Cj,,km、Cjkm,,、Cjkm,,分别表示图像在水平、垂直、对角方向上的高频分量。2.对小波分解得到每一级高频系数的各个方向的系数作对应方向的方向检测和53 InSAR干涉纹图噪声抑制方法研究中值滤波:(1)检测边缘方向。对不同级别各个方向的小波高频系数,以每个高频系数为中心,用大小为w×w的窗口,其中,w=5+2×(k-1),k为小波分解级数,取kkk值1~3,对各级各个方向的高频系数,分别检测此方向系数在各个位置处的各个方向强度变化程度,由此判断此系数是否为对应方向边缘。对位置(i,j)处系数00p(,)ij,首先令w=int(w/2)-1,int函数表示取整数部分,计算高频系数各个方向00k的平均强度。水平方向:jw0+Rp(1,1)=(∑(,))/(i0jwk−2)jjw=−0jw0+Rp(1,2)=(∑(i0−−1,))/(jwk2)(5.17)jjw=−0jw0+Rp(1,3)=(∑(i0+−1,))/(jwk2)jjw=−0垂直方向:iw0+Rp(2,1)=(∑(,ij0))/(wk−2)iiw=−0iw0+Rp(2,2)=(∑(,ij0−−1))/(wk2)(5.18)iiw=−0iw0+Rp(2,3)=(∑(,ij0+−1))/(wk2)iiw=−0045对角线方向:iwjw00−+,Rp(3,1)=(∑(,))/(ijwk−2)iiwjjw=+00,=−iwjw00−++,(1)Rp(3,2)=(∑(,))/(ijwk−1)(5.19)iiw=++00(1),jjw=−iwjw00−++(1),Rp(3,3)=(∑(,))/(ijwk−1)iiwjjw=+00,=−+(1)0135对角线方向:iwjw00++,Rp(4,1)=(∑(,))/(ijwk−2)iiwjjw=−00,=−iwjw00+++(1),Rp(4,2)=(∑(,))/(ijwk−1)(5.20)iiwjjw=−00,=−+(1)iwjw00+++,(1)Rp(4,3)=(∑(,))/(ijwk−1)iiw=−+00(1),jjw=−54 硕士学位论文然后计算此位置各个方向系数强度变化的较大值。水平方向:RRRRR=−max((1,1)(1,2),(1,1)−(1,3))(5.21)1垂直方向:RRRRR=−max((2,1)(2,2),(2,1)−(2,3))(5.22)2045对角线方向:RRRRR=−max((3,1)(3,2),(3,1)−(3,3))(5.23)30135对角线方向:RRRRR=−max((4,1)(4,2),(4,1)−(4,3))(5.24)4令RR=max(,R,R,R)1234,为各方向最大强度变化。再对某个方向系数判断R是否为此方向的R,若R=R,此系数为表示此方向的边缘,否则,为此方向高频rr系数的非边缘。(2)如果第l级、r方向的边缘位置(i,j)处的系数为此方向系数,同时进一步00根据此方向R的值,判断方向的取向,例如对水平方向的边缘判断是水平方向上r边缘或下边缘等等,最后根据下面公式得到滤波后此位置处对应的系数pi(,)j。lr,00pijm(,)=edianP(DMg)D(5.25)lr,00l,rl,r式中,median为中值运算,PDl,r为以边缘处(,)ij00为中心,窗口大小wwkk×、第l级、r方向的高频系数数据,MD第l级、r方向的模板,边缘的不同方向取相应级l,r别窗口大小的不同方向的模板,例如第一级窗口的方向模板如图5.15所示。(3)对第l级、r方向的非边缘位置(,)ij处系数,根据公式(5.26)得到滤波后此00位置处对应的系数值p(,)ij。lr,00pijm(,)=edianP(D)(5.26)lr,00l,r式中,PDl,r为以(,)ij00为中心,以窗口大小(wk+2)×(wk+2)、第l级、r方向的高频系数。3.由分解后的低频部分系数和方向检测中值滤波处理后的高频系数重构InSAR干涉图。55 InSAR干涉纹图噪声抑制方法研究0000011111001111100000000111110011111000111110000000111110001111100000001111100011111000000011111000a)水平方向1b)水平方向2c)垂直方向1d)垂直方向20000111111100001111100011111101100001111001111110011100001110111111000111100001111111100001111100001。。。。e)45方向1f)45方向2g)135方向1h)135方向2图5.15滤波窗口模板5.5.2实验结果与分析基于小波变换和方向中值滤波方法的实验结果如图5.16和5.17所示。这里选择可以应用于信号和图像重构的正交小波函数bior5.5,用ATrous算法实现小波分解,小波分解的层数N为3,窗口大小选取为5×5。a)喀什地区原始干涉纹图b)小波变换和方向中值滤波后的图像图5.16喀什地区小波变换和方向中值滤波效果图56 硕士学位论文a)乌鲁木齐地区原始干涉纹图b)小波变换和方向中值滤波后的图像图5.17乌鲁木齐地区小波变换和方向中值滤波效果图实验前后干涉纹图中的残差点比较如表5.4所示。表5.4小波和方向中值滤波残差点比较表噪声抑制方法残差点数喀什地区原始干涉纹图(无滤波)54428小波变换和方向中值滤波(5×5)1054乌鲁木齐地区原始干涉纹图(无滤波)79684小波变换和方向中值滤波(5×5)2036实验证明,在小波变换的基础上,通过对各个方向的高频系数中对应方向的边缘位置平滑后做中值滤波处理来实现干涉图的滤波,滤波后的干涉纹图视觉特征良好,残余点大大减少,很好地保持了边缘特征信息。57 InSAR干涉纹图噪声抑制方法研究结论干涉合成孔径雷达(InSAR)三维成像,是合成孔径雷达遥感技术的新发展,是一种全新的对地观测与感知技术。它通过两副天线同时观测或通过一副天线两次平行观测,获取地面同一景观的复图像对(包括强度信息和相位信息),根据地面各点在两幅复图像中的相位差,得出各点在两次成像中微波的路程差,计算出所摄图像地区的地形、地貌以及表面的微小变化。这一技术可以广泛应用于军事、地质、大地测量、地震、地质灾害监测等领域。本论文的研究工作在国家自然科学基金(60375001)和教育部重大项目培育基金(706043)的资助下,根据InSAR三维成像技术的特点,从InSAR系统的处理过程出发,分别对InSAR系统的工作原理、InSAR干涉纹图噪声特性、InSAR干涉相位纹图噪声抑制方法进行了深入研究。主要工作如下:1.总结了InSAR三维成像技术的发展历史和现状,详细研究了InSAR三维成像技术的原理。介绍了InSAR技术的三种干涉模式和它们的特点,以重复轨道干涉模式为例,阐述了InSAR三维成像技术的原理,给出了相应的数据处理流程,并对各个处理步骤作了简要说明,同时分析了影响InSAR系统测量精度的主要因素。2.对InSAR干涉纹图的噪声来源、特性和相位统计特征进行了深入地分析。3.对传统的干涉纹图滤波方法进行了分析和总结,以实验结果为基础,全面分析了这些算法的优缺点。4.将小波分析引入干涉纹图噪声抑制当中,提出了一种基于小波变换和加权圆周期中值的干涉纹图滤波方法。在综合考虑高频系数方向特性和统计特性的基础上,研究了一种基于小波变换和方向中值的干涉纹图滤波方法,该方法在小波变换的基础上,分别对各个方向的高频系数部分检测其是否为对应方向的边缘,对不同的边缘方向用不同的模板平滑后再进行中值滤波处理。5.利用德国欧洲空间局(ESA)提供的两对ERS-1/2串行SAR数据进行了实验,实验结果验证了本文所提出的处理算法的可靠性及有效性。本论文的研究工作为以后的数据处理和工作奠定了基础,但是还不够全面,许多方面还需要进一步的研究,作者认为有以下几个方面的工作可以进一步开展:(1)研究基于神经网络和模糊神经网络的干涉纹图噪声抑制方法;(2)研究基于多尺度形态学操作的干涉纹图噪声抑制方法;(3)研究多滤波方法数据融合的干涉纹图噪声抑制方法。58 硕士学位论文参考文献[1]WMBrownm.Syntheticapertureradar.IEEETrans.Aerosp.Elec.Syst,1967,3(2):217-229[2]CAWiley.Syntheticapertureradar.IEEETrans.Aerosp.Elec.Syst.,1985,21(3):440-443[3]RBirk,WCamus,EValenti,etal.Syntheticapertureradarimagingsystem.IEEEAESSystemMagazine,1995,10(11):15-23[4]ZebkerHA.StudyingtheEarthwithinterferometricradar.ComputinginScience&Engineering,2000,2(3):52-58[5]LiFK,RMGoldstein.Studiesofmultibaselinespaceborneinterfero-metricsyntheticapertureradars.GeoscienceandRemoteSensing,IEEETransactionson,1990,28(1):88-94[6]LCGraham.Syntheticinterferometricradarfortopographicmapping.Proc.ofIEEE,1974,62(6):763-768[7]RosenPA,etal.Syntheticapertureradarinterferometry.ProceedingsoftheIEEE,2000,88(3):333-342[8]RodriguezE,JMMartin.Theoryanddesignofinterferometricsyntheticapertureradars.RadarandSignalProcessing,IEEProceedingsF,1992.139(2):147-155[9]舒宁.雷达影象干涉测量原理.武汉:武汉大学出版社,2003[10]HenriMaitre编,孙洪等译.合成孔径雷达图像处理.北京:电子工业出版社,2005[11]HAZebker,RMGoldstein.TopographicMappingFromInterfero-metricSyntheticApertureRadarObservations.JournalofGeophysicalResearch,1986,91(B5):4993-4999[12]MocciaA,Verella.Atetheredinterferometricsuntheticsyntheticapertureradarforatopographicmission.IEEETrans.Geosci.andRemotesensing.1992,30(1):103-109[13]AEERogersandRPIngalls.Venus:mappingthesurfacereflectivitybyradarinterferometry.Science,1969,165(8):797-799[14]GoldsteinRM,BarnettTP,ZebkerHA.Remotesensingofoceancurrent.Science,1989,24(6):1282-128559 InSAR干涉纹图噪声抑制方法研究[15]GuarnieriAM,PratiC.SARinterferometry:A"QuickandDirty"coherenceestimatorfordatabrowsing.IEEETransactiononGeosciencyandRemoteSensing.1997,35(3):660-669[16]FKLi,RMGoldstein.StudiesofMultibaselineSpaceborneInterferometricSyntheticApertureRadars.IEEETrans.GeoscienceandRemoteSensing,1990,1(28):88-97[17]FabioGatelli,MontiGuarnieri,etal.TheWavenumberShiftinSARInterferometry.IEEETrans.Geosci,RemoteSensing,1994,32(4):855-865[18]ZebkerHA,VillasenorJ.DecorrelationinInterferometricradarechos.IEEETrans.GeoscienceandRemoteSensing,1992,30(5):950-959[19]GrahamLC.Syntheticradarfortopographicmapping.Proc.IEEE,1974,7(62):763-768[20]GabrielAK,GoldsteinRM,Crossedorbitinterferometry:theoryandexperimentalresultsfromSIR-B.lnternationalJournalofRemoteSensing,1988,9(8):852-857[21]LeeJS,PapathanssiouKP,etal.ANewTechniqueforNoiseFilteringofSARInterfermetricPhaseImages.IEEETransonGeoscienceandRemoteSensing,1998,36(5):1456-1465[22]RBamler,DJust,Phasestatisticsininterferogramswithapplicationstosyntheticapertureradar,AppliedOptics,1994,33(20):4361-4368[23]BaranL,etal.AmodificationtotheGoldsteinradarinterferogramfilter.GeoscienceandRemoteSensing,IEEETransactionson,2003,41(9):2114-2120[24]Jong-SenL,etal.IntensityandphasestatisticsofmultilookpolarimetricandinterferometricSARimagery.GeoscienceandRemoteSensingIEEETransactionson,1994,32(5):1017-1026[25]RMGoldstein,HAZebker,CLWerner.SatelliteradarinterferometryTwo-dimensionalphaseunwrapping.RadioScience,1998,262(6):1525-1530[26]王志勇,张继贤,黄国满.InSAR干涉条纹图去噪方法的研究.测绘科学,2004,12(29):31-33[27]A.K.Gabriel,R.M.GoldsteinandH.A.Zebker.Mappingsmallelevationchangesoverlargeareas:Differentialradarinterferometry.J.Geophys.Res,1989,94(B7):9183-9191[28]郭华东,徐冠华.星载雷达应用研究.北京:中国科学技术出版社,1996[29]王超,张红,刘智.星载合成孔径雷达干涉测量.北京:北京科学出版社,2002[30]王超.利用航天飞机成象雷达干涉数据提取数字高程模型.遥感学报,1997,60 硕士学位论文1(1):46-49[31]潘习哲.星载SAR图象处理.北京:科学出版社,1996[32]郭华东.航天多波段全极化干涉雷达的地物探测.遥感学报,1997,1(1):32-39[33]EinederM,SSuchandt.RecoveringradarshadowtoimproveinterferometricphaseunwrappingandDEMreconstruction.GeoscienceandRemoteSensing.IEEETransactionson,2003,41(12):2959-2965[34]TrouveE,JMNicolas,HMaitre.Improvingphaseunwrappingtechniquesbytheuseoflocalfrequencyestimates.GeoscienceandRemoteSensing.IEEETransactionson,1998.36(6):1963-1972.[35]MBao,CBruning,WAlpers.Simulationofoceanwavesimagingbyanalong-trackinterferometricsyntheticapertureradar.IEEETrans.Geosci.RemoteSensing,1997,35(3):618-631[36]CIchoku,AKaenieli,YArkin,etal.ExploringtheutilitypotentialofSARinterferometriccoherenceimages.Int.J.Remotesensing,1998,19(6):1147-1160[37]RNTreuhaft,SNMadsen,MMoghaddam,etal.Vegetationcharacteristicsandunderlyingtopographyfrominterferometricradar.RadioScience,1997,31(6):1449-1485[38]LanariR,FornaroG,etal.GenerationofDigitalElevationModelsbyUsingSIR-C/X-SARMultifrequencyTwo-passInterferometry:theEtnaCaseStudy.IEEETrans.onGRS,1996,34(5):1097-1114.[39]徐华平,周荫清,陈杰等。干涉SAR中相位图的噪声抑制.北京航空航天大学学报,2001,2(27):16-19[40]JOHagberg,LMHUlander,JAskne.Repeat-passSARinterferometryoverforestedterrain.IEEETrans.Geosci.RemoteSensing,1995,33(2):331-340[41]崔锦泰著,程正兴译.小波分析理论.西安:西安交通大学出版社,1998.[42]飞思科技产品研发中心编著.小波分析理论与MATLAB7实现.北京:电子工业出版社,2005.[43]杨福生.小波变换的工程分析与应用.北京:科学出版社,1999.[44]程正兴.小波分析算法与应用.西安:西安交通大学出版社,1998.[45]胡昌华,张军波,夏军等.基于MATLAB的系统分析与设计——小波分析.西安:西安电子科技大学出版社,1999.[46]孟丹.InSAR干涉相位图的噪声滤波:[吉林大学硕士学位论文],2006:18-41[47]段克清,向家彬,汪枫.InSAR相位条纹图的加权圆周期中值滤波算法.空61 InSAR干涉纹图噪声抑制方法研究军雷达学院学报,2005,3(1):4-6[48]张恒,雷志辉,丁晓华.一种改进的中值滤波方法.中国图像图形学报,2004,9(4):408-411[49]林卉,赵长圣,杜培军等.InSAR干涉图滤波方法研究.测绘学报,2005,5(34):113-117[50]穆冬,朱兆达,张焕春.干涉SAR相位条纹的鲁棒加权圆周期滤波.数据采集与处理,2001,9(16):299-303[51]于晶涛,陈鹰.一种新的InSAR干涉条纹图滤波方法.测绘学报,2004,5(33):121-126[52]孙龙,胡茂林,张长耀.一种新的InSAR干涉纹图的滤波方法.合肥工业大学学报(自然科学版),2005,1(28):79-83[53]单世铎.InSAR平地效应去除及相位噪声抑制方法研究:[解放军信息工程大学硕士学位论文],2005:11-17[54]韩春明等.一种改进的SAR图像斑点噪声滤波方法.遥感学报,2004,8(2):121-12762 硕士学位论文附录A攻读学位论文期间所发表的学术论文目录[1]刘实,毛建旭.基于小波变换和加权圆周期中值算法的InSAR干涉纹图除噪[J].中国仪器仪表,2007,03:56-59[2]刘实,毛建旭.基于小波变换和圆周期均值算法的InSAR干涉纹图噪声抑制[J].测控技术,2007,第10期,已录用,待发表.63 InSAR干涉纹图噪声抑制方法研究致谢值此论文完成之际,首先我要对我的导师毛建旭老师致以衷心的感谢。本文的撰写是在毛老师的悉心关怀和精心指导下完成的。感谢导师三年来在学习和工作上对我的启发和教益,在生活上对我的关心和照顾,令我终生难忘。毛老师渊博的学识、忘我的工作热情深深感染着我。本文从选题到理论研究,直到论文的撰写,每一阶段无不凝聚着导师的无私奉献和辛勤劳动。同时,毛老师性格随和,待人宽容有礼,对待学术严谨求实,这些都是我在今后工作学习中的榜样。衷心感谢毛老师对我的培养。在电气与信息工程学院的三年学习与生活中,衷心感谢王耀南老师、戴瑜兴老师、孙炜老师、王绍源老师、彭楚武老师、李树涛老师、王玲老师等的无私关心和帮助。在遥感实验室的学习期间,实验室的余洪山博士、彭金柱博士、何儒云博士、余浩硕士、以及同届的张辉、王威、刘玲、李红英、蔡俊伟等给了我很多帮助,尤其学术上的启发让我获益匪浅。同时要感谢实验室的周博文、张志国、蔡玉连和师弟刘洋、岑小林等,让我这段校园生活留下了很多美好的回忆。最后,我要感谢我的父母、家人给予我生活和学业上的支持和帮助,你们是我前进道路上最坚强的后盾!刘实2007年4月64