《insar相位解缠算法分析与其软件开发》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
西南交通大学硕士研究生学位论文第1I页相位数据测试了InSAR相位解缠系统的执行效率和稳定性,同时评定了系统提供的八种解缠算法的精度,得出各相位解缠算法的特性和优缺点如下:(1)基于离散余弦变换(DCT)的最小二乘算法所需内存最少,执行速度也最快,但它与无权重多重网格算法一样,解缠结果出现与初始缠绕相位不一致并且更加平滑的机率较高。因此,无权重最小二乘解缠算法在实际应用中适用性较差,一般不建议使用。(2)Goldstein“枝切"算法执行速度较快,占用内存较少,对信噪比高的相位解缠精度较高,在实际应用中常作为首选。但对于干涉图质量差异较大的区域,解缠结果容易出现不连续或许多孤立区域。然而利用系统提供的手动修改“枝切"线和框选复解缠操作可改善解缠效果,使实验中的“mountains”解缠相位标准偏差从0.099rad减d,N0.047tad,相应的DEM差异的标准偏差从3.409m降低到1.698m,精度得到了很大提高,可以和效果最好的最小F范数算法相媲美,但其运算速度却比最小矽范数算法快很多。(3)当有好的质量图可利用时,可选用质量图路径引导算法、“掩膜枝切’’算法以及预处理共轭梯度算法,若是前几种算法均解缠失败,最后才考虑解缠效果较好的Flynn最小不连续算法和最小r范数算法。尤其是最小F范数算法,解缠相位与SARseape最小费用流算法的解缠相位差异均值为3.119rad,标准偏差仅为0.041rad,相应的DEM差异均值为0.087m,标准偏差为1.439m。但是,Flynn最小不连续算法和最小∥范数算法所消耗的内存也较多,执行速度也较慢,其中最小r范数算法的速度最慢。(4)总的来讲,最小二乘算法解缠结果较平滑,对噪声区域也能解缠,但同时也将误差传播到高质量区域,造成系统性偏差;路径跟踪算法将噪声区或相位梯度变化大的区域进行隔离,让这部分相位不参与解缠,从而避免了误差的传递。因此,在相位信噪比高和梯度变化小的区域可选用路径跟踪相位解缠算法,信噪比低或梯度变化大的地方则可使用最小二乘相位解缠算法。关键词:合成孔径雷达干涉测量;相位解缠;路径跟踪;最小范数;系统开发 西南交通大学硕士研究生学位论文第1II页AbstractSyntheticApertureRadarIntefferometry他S趾◇andDifferentialInSAR①-InsAR)aremicrowaveremotesensingtechnologiesdevelopedsincethe1960s.Comparedwiththetraditionalremotesensingtechnology,InSARhasmanytechnologicaladvantages,suchaslargecoverage,workingatalltimeandunderall-weatherconditions.Ithasbeen埘delyusedintheextractionofdigitalelevationmodel(DEM)andthesurfacedeformationmeasurements.PhaseunwrappingisoneofthekeystepsandtechnicaldifficultiesinthedamprocessingofInSARandD-InSAR.ItsaccuracywilldirectlyaffecttheprecisionoftheDEMextractedbyInSARandthegrounddeformationmeasuredbyD-InSAR.Inthepast30years,scholarsathomeandabroadhaveproposedmanyrepresentativealgorithmsofphaseunwrapping.Accordingtothestrategiesofunwrappingprocedures,theyCanbedividedintotwocategories:thepath-followingmethodandthemini/nunlnormmethod.Inspecificapplications,phaseunwrappingismainlydesigned嬲anembeddedmoduleincludedintheInSARandD-InSARdataprocessingsol,wares,andthereisnosuchprofessionalsoRwarewhichisspecificallydevelopedforphaseunwrapping.Moreoverthecurrentsoft:waresonlyCOVeraroundonetotwokindsofunwrappingmethods,whichresultsinpoorselectivityandadaptability.Furthermore,currentresearchesonphaseunwrappingoftenfocusonimprovingtheunwrappingalgorithms,andthereislackofcomprehensivecomparationanalysisandevaluationofdifferentalgorithms.Tllisarticleemphasisonstatingtheprocessingproceduresandimplementationmethodsoffourpath-followingmethods一(Goldstcinbranchcutalgorithm,quality-guidedpath-followingalgorithm,maskcutalgorithmandFlynn’Sminilnunldiscontinuityalgorithm)andfourminilnulnnormmethods(discretecosinetransformalgorithm,unweightedmultigridalgorithm,preconditionedconjugategradientalgorithmandminin】umLe-normalgorithm)basedonintroductionofthebasicprinciplesofphaseunwrapping.OntheWindowsoperatingsystem,combiningwiththeobject-orientedC撑languageandArcEnginewhichistheembeddedcomponentlibraryofGIS,aInSARphaseunwrappingsystemwhichisdedicatedtophaseunwrappingWasdeveloped.ToimprovetheSuccessrateandreliabilityofphaseunwrapping,theman-machineinteractiveoperationandtheevaluationstrategywereintroducedintheprocessofphaseunwrapping,i.e.,addmodulestothesystemwhichincludeunwrappingevaluationmodule,DEMgenerationmodule,moduleofmanuallymodifytheunwrappingpath(i.e.,thebranchlines) 西南交通大学硕士研究生学位论文第Ⅳ页andboxselectedunwrappingmodule.Inaddition,toimprovetheeffectofunwrapping,thesystemalsoprovidesthesub-moduleofimageprocessing.ExceptGAMMAsoftware,theInSARphaseunwrappingsystemdefinesthesanledatastructureasotherInSARandD-InSARprocessingsoRwareswhicharecommonlyused,andthusitallowsdatasharingwithothersoRwares.What’smore,thesystemprovideseightkindsofphaseunwrappingmethods,SOuserscallselecttheappropriateunwrappingmethodaccordingtodifferenttypesofdata.Withtheaidoftheunwrappingauxiliarymodulestoimprovetheunwrappingresults,theprecisionoftheeventuallyacquiredDEMandgrounddeformationscallbepromoted.Finally,boththesimulatedwrappedphasedataandtheinterferometricphasedataofsoutheastLincolncountyofUnitedStateswereusedtotesttheexecutionefficiencyandstabilityoftheInSARphaseunwrappingsystem,andtoevaluatetheaceuracyoftheeightalgorithmsinthesystem.11圮characteristics,advantagesanddisadvantagesofthesealgorithmsareasfollows:First,least-squaresalgorithmbasedonthediscretecosinetransform(DCT)reqlliI订theleastamountofmemo巧and晰tllthefastestexecutionspeed.Butthereishigherprobabilityofresultingininconsistenciesbetweentheresultsandtheinitialphaseandtheresultsoftenappearssmoother,whichissimilar.withtheunweightedmultigridalgorithm.Therefore,unweightedleast-squaresmethodhaspoorapplicabilityinpractice,SOitisgenerallynotrecommendedtobeused.Second,Goldsteinbranchcutalgorithmhashighexecutionspeedandlowmemo巧needs.Theprecisionofunwrappingresultishigherwhenthephasehashi曲signal-to-noiseratio,SOitisselectedinthefirstplaceinthepracticalapplication.Butfortheareaswherethequalityofinterferencemapsvaryremarkably,theresultsofunwrappingarepronetobediscontinuousormanyisolatedareasmaybeinduced.However,afterusingofthemanuallymodifiedbranch··cutlinesandboxselectedre·-unwrappingoperationwhichprovidedbythesystem,theeffectofunwrappingCanbeimproved,theymadethestandarddeviationoftheunwrappingphaseofthe”mountains”intheexperimentdecreasesfromO.099radto0.047tadandthestandarddeviationofthecorrespondingDEMdecreasesfrom3.409mto1.698m.Theaccuracyhasbeengreatlyimproved,andtheeffectcanbecomparablewiththeminimll/nrnormalgorithm,butitsoperationspeedismuchfasterthantheminimumtnormalgorithm。删,whengoodqualitymapisavmlable,wecanchoosequality-guidedpath—following 西南交通大学硕士研究生学位论文第V页algorithm,maskcutalgorithmandpreconditioningconjugategradientalgorithm.Ifthesealgorithmsareallfailure,thentheFlynn’SminilnunldiscontinuityalgorithmandmlnllnaInLe-Normalgorithmcanbeconsidered.Especially,themeanofdifferencesbetweentheunwrappedphasefromminJlnuinif-NormalgorithmandminimumcostflowalgorithminSARscapeis3.119tad,andthestandarddeviationisonly0.041rad,themeanofthedifferencesbetweenthecorrespondingDEMis0.087mandthestandarddeviationis1.439m.However,Flyun’S蛐UlXIdiscontinuityalgorithmandminlnlunlLe-Normalgorithmconsumemorememories,andtheexecutionspeedisslowerwithminilnulTlLP-Normalgorithmistheslowestone.Forth,ingeneral,theunwrappingresultsofleast-squaresmethodsaresmoother,itcanunwrapthephaseinthenoisyareas,butatthesametime,theerrorswillpropagatetothehighq砌时areas,resultinginsystemicbias.Path-followingmethodswillisolatethenoiseareasandtheregionswherethechangeofphasegradientislargeanditlimitsthephaseinsuchareastoparticipateintheunwrapping,thusCanavoidthetransmissionoferrors.Asaresult,inareaswherethesignal-to—noiseratioishighandgradientchangeissmall,wecanchoosePath-followingmethods,whileinareaswithlowsignal-to-noiseratioorhighphasegradientchangeratetheleast-squaresmethodsseemmoreappropriate.Keywords:InSAR;PhaseUnwrapping;PathFollowing;MininluRlNorm;SystemDevelopment 西南交通大学硕士研究生学位论文第VI页目录第1章绪论⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯11.1研究背景及意义⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯11.1.1研究背景⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。l1.1.2研究意义⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯31.2InSAR相位解缠的国内外研究现状⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.41.2.1相位解缠算法的研究现状⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。41.2.2相位解缠软件研发及应用现状⋯⋯⋯⋯⋯⋯⋯⋯j⋯⋯⋯⋯⋯⋯⋯⋯。51.3论文研究内容⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.61.4论文组织结构⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯7第2章相位解缠原理及关键因素分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..92.1引言⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯92.2基于路径跟踪的相位解缠基本原理⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.92.2-1留数定理⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.102.2.2留数点探测⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.1l2.3基于最小范数的相位解缠基本原理⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..122.4关键因素分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.142.4.1质量图⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..142.4.2掩膜⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..152.4.3噪声滤波⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯152.5解缠评价标准⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯l62.5.1误差图⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..162.5.2不连续图⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.172.5.3反缠绕相位图⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..172.6InSAR解缠相位的应用—.DEM的生成⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。172.6.1生成高程图⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯172.6.2地理编码⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯l82.7小结⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.19第3章相位解缠算法的实现⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..20 西南交通大学硕士研究生学位论文第VII页3.1引言⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.203.2基于路径跟踪的相位解缠算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..203.2.1Goldstein“枝切’’算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.203.2.2质量图路径引导算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.233.2.3“掩膜枝切"算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.243.2.4Flyrm最小不连续算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯263.3基于最小范数的相位解缠算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯293.3.1无权重最小二乘算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.293.3.2加权最小二乘算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.333.3.3最小r范数算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯333.4小结⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.35第4章InSAR相位解缠程序设计与开发⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯364.1引言⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯364.2系统开发平台⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯364.2.1操作系统⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯364.2.2开发工具⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.364.3系统的开发与实现⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。374.3.1系统主界面⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯374.3.2数据输入界面⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.384.3.3数据处理界面⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.384.3.4框选解缠操作⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯:⋯⋯⋯⋯⋯⋯⋯⋯444.3.5手动修改“枝切"线操作⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯474.4,J、结⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.49第5章基于模拟数据的解缠实验与结果分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯505.1引言⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..505.2实验数据⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯505.2.1模拟InSAR数据⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯505.2.2“剪切面"数据⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯515.2.3受高斯噪声干扰的“剪切面’’数据⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.51 西南交通大学硕士研究生学位论文第VIII页5.3相位解缠实验与结果分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。525.3.1Goldstein“枝切”算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..525.3.2质量图路径引导算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..555.3.3“掩膜枝切”算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.575.3.4Flyrm最小不连续算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯595.3.5无权重最小二乘算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯625.3.6预处理共轭梯度算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯635.3.7最小r范数算法⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯665.4小结⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.68第6章相位解缠算法对比分析与精度评定⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯706.1引言⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.706.2研究区域和实验方案⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。706.2.1研究区域⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。706.2.2实验方案⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯706.3相位解缠及DEM生成⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯706.3.1自适应滤波⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯..706.3.2相位解缠⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。726.3.3生成高程图与地理编码⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯736.4解缠算法精度评定与对比分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯736.4.1精度评定⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯736.4.2对比分析⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯。766.5小结⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.79总结与展望⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯80致谢⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯83参考文献⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯.84 西南交通大学硕士研究生学位论文第1页第1章绪论1.1研究背景及意义1.1.1研究背景合成孔径雷达(SyntheticApertureRadar,SAR)作为一种高分辨率的二维成像雷达,与传统的可见光、红外遥感成像传感器相比,具有显著的优越性。SAR属于主动式微波传感器,其成像不依赖于太阳辐射能量,通过自身对目标发射脉冲波束并接收其后向散射信号进行成像,且微波能够穿透云雾、烟尘、雪地以及被沙漠或植被覆盖的地质地物等【¨】。这使得SAR具有全天候、全天时的对地观测能力,可弥补光学传感器在时间与空间上受限制造成的成像“盲区”,有效获取地面的高质量高分辨率影像【41。合成孔径雷达干涉技术(SyntheticApertureRadarInterferometry,:nSAR)是以合成孔径雷达复数数据提取的相位信息为基础来获取地表三维信息和形变信息的一种技术[51。InSAR通过两副天线同时观测(单轨模式),或两次近平行的观测(重复轨道模式,如图1.1所示),获取地面同一景观的复数影像对K堋。将两幅SAR复数影像中对应像素的复数值共轭相乘(即相位值相减),可得到干涉相位图。该干涉相位的取值在(一乃7/"】之间,.万至万的周期性变化呈现为一个干涉条纹【9,.0l。由于SAR复数影像中的相位信息包含了天线和地表距离信息,所以可以对覆盖同一地区的两幅SAR影像的相位进行干涉处理,并利用其雷达平台的姿态参数重建地表三维信息,即数字高程模型(DEM),可达到米级的精度【11.131。InSAR的数据处理流程如图1-2所示。图1-1InSAR的成像几何关系 西南交通大学硕士研究生学位论文第2页图1—2InSAR获取数字高程模型(DEM)的流程图为了得到地表高程,必须为干涉相位图中每个缠绕的相位加上或减去2万的整数倍,即相位解缠。相位解缠是合成孔径雷达干涉数据处理中的重要步骤,它的准确性将影响到最终提取的数字高程图的精度。从图1.1所示的雷达成像几何关系可以推出高程的计算公式:让sin-1(篙)川(1-2)h=H一蜀cos0(1.3)其中,除绝对相位矽外,其余参数均可从SAR头文件中的轨道姿态数据中推算出来,而矽只能以相位解缠的方式计算出来【141。合成孔径雷达差分干涉技术(DifferentialInSA艮D.InSAR)是对InSAR技术应用的一个扩展,主要用来监测地球表面的微小形变,可达厘米、甚至毫米级的精度‘7,81。D.InSAR广泛应用于火山运动、地震形变、冰川漂移、地面沉降及山体滑坡等的监测【15一18】oD.InSAR同样对不同时间所获取的两幅SAR影像(成像期间地表发生了形变)做干涉处理,得到的干涉相位不仅包含了地形的高程信息,而且还包含了两次成像期间的地表形变信息‘191。从干涉相位中去除反映高程的相位,即可得到反映成像期间地表形变的相位。但该相位的值域仍在(听,万】之间,因此也必须经过相位解缠求出其真实相 西南交通大学硕士研究生学位论文第3页位,才能计算出雷达斜距向上的地表形变量。D.InSAR数据处理流程如图1.3所示。{幽|一1.1.2研究意义豳图1.3D.InSAR获取地表形变量的流程图(两轨法)l渗:无论是应用InSAR提取地表高程还是应用D.InSAR获取地表形变的数据处理过程中,相位解缠都是最重要、最关键的步骤之一,相位解缠的准确性直接影响到生成的高程和形变测量的精度[20之61。然而,各种干扰因素不可避免的影响着SAR干涉图数据,使其出现不连续以及干涉条纹紊乱、不清楚等现象,使得相位解缠的进行变得十分艰难[27。11。值得指出的是,虽然目前关于InSAR和D.InSAR数据处理的软件均提供了相位解缠算法,但算法均较单一,而且解缠过程封装性较强、可视化程度较低、交互式操作较少。此外,相位解缠算法种类虽多,但普适性太差,没有一种通用的算法,并且不同的解缠算法在解缠效率、精度等方面存在不同,对于解缠相位精度的评估不够全面、充分,仍没有一个精度评估的标准。因此,为了弥补当前InSAR相位解缠研究及应用中所呈现的不足,有必要对相位解缠的原理和算法进行全面、深入的分析,并对相应的算法进行实现,形成一个具有可视化、人机交互式操作特点并且解缠算法较全面的相位解缠系统,同时对各种解缠算法进行评定与精度验证。 西南交通大学硕士研究生学位论文第4页1.2InSAR相位解缠的国内外研究现状1.2.1相位解缠算法的研究现状20世纪60年代末至70年代末,以采用积分法为理论基础的一维相位解缠发展到二维。随后,由于InSAR等二维图像处理的需要,二维相位解缠技术迅速发展。近30年来,为了寻求高质量的相位解缠结果,国内外学者做了大量研究并提出了许多具有代表性的相位解缠算法,按解缠方式可分为两类:(1)基于路径跟踪(Path-Following)的相位解缠算法,利用相邻相位在空间上的相关性进行解缠;(2)基于最小范数(Minimum-Norm)的相位解缠算法,将相位解缠问题转化为数学上的最小范数问题来解决。基于路径跟踪相位解缠算法最早由Goldstein等人于1988年提出,后来成为经典相位解缠算法一“枝切"(Branch-Cut)'法1321。该算法通过定义、识别留数点,依据最近邻原则连接留数点形成“枝切”线来隔离噪声、防止局部误差传播到整个图像中,然后对相邻像元的缠绕相位梯度进行积分,从而实现相位的解缠【33‘34]。对于数据质量好的局部相位的解缠精度较高,并且解缠速度快;但当“枝切’’线放置不合理时将会导致误差的传播。这是由于“枝切”法只利用了“留数点"信息,所以通常会有“枝切"将大部分区域隔离起来形成闭合回路以及“枝切"安置错误的现象。在Goldstein的基础上,Huntley与Cusaek对枝切线的连接进行了改进,使其对噪声更具免疫性[35,36]。1991年,Bone提出利用质量评价来指导相位解缠,对相位数据求二次偏导,设置阈值作为象元的质量标准,仅对相位值的二次偏导低于阈值的相位进行解缠[371。而Flyrm在1996年提出了用质量图来指导连接枝切线的新思想,这种算法被称为掩膜枝切法(M嬲k-cut)【3引。1997年,Flyrm等人又提出了Flynn最小不连续的相位解缠算法,该算法通过使用网络图来选择最优积分路径,从而使干涉相位中不连续的长度达到最d4391。1999年,WeiXu等人提出的区域增长算法主要是使用外部辅助信息——质量图来引导积分路径的选择【柏】。国内学者岑等在2008年提出了将质量图与留数点相结合的相位解缠算法【4l】。同年,谢等人基于枝切截断又提出了一种可自动避开不连续相位的算法[421。最小二乘相位解缠算法最早由Frid等人于1977年提出,同时分析了最dx--乘算法解缠精度的问题,并提出一种基于噪声传播最小准则下的线性优化算法【43】。但是这 西南交通大学硕士研究生学位论文第5页种算法没有权重,因此相位噪声对解缠结果的影响较大。1994年,Ghiglia与Romero提出了一种加权最小二乘解缠算法——预处理共轭梯度法(PreconditionedConjugateGradientMethod,PCG),其基本思想是先用无权重最小二乘算法(DCT)作为预处理条件改善加权正态方程系数矩阵的条件数,加速收敛;再采用共轭梯度法迭代获得绝对相位Ⅲ】。后来,Pritt等对此方法做了改进,将快速傅里叶变换ffastFourierTransform,FFT)引入不加权最小二乘相位解缠中,解决了泊松方程的边界问趔451。1996年Ghiglia与Romero针对相位解缠问题的一般r范数解提出了最小r范数算法,该算法可以不依赖于质量图,根据干涉相位自身便可生成数据相关的权重【46】。1998年,Costantini提出了基于网络规划的相位解缠算法,将最小化问题转为求解最小费用流的网络优化问题,很好的解决了相位解缠效率和解缠准确性不兼顾的问题【47】。2007年,陈提出了在粗细大小不同的网格上进行迭代计算的多重网格算法【48】。陈等在2008年又通过使用相位导数变化图来定义权重,提出了基于小波变换的算法【49】。每种算法都有其自身的优缺点,适用于特定条件的数据,普适性都不很好,因此应根据实际情况来选择合适的解缠算法。1.2.2相位解缠软件研发及应用现状在过去20几年间,合成孔径雷达干涉测量技术从数据获取、技术研发到应用逐步走向成熟。SAR数据的处理从早期需要过多的人工参与到现在基本实现了自动化【5】,极大的提高了InSAR数据处理的效率,同时也伴随出现了许多InSAR数据处理的商业软件和科研软件。但是,到目前为止相位解缠主要是作为一个嵌入模块被包含在InSAR和D.InSAR数据处理软件中,专门针对相位解缠的专业性软件极少。相应的InSAR和D.InSAR数据处理软件主要有【50。53】:(1)基于Unix或Linux或WindowsCygwin操作系统的GAMMA、DORIS和ROI—PAC;(2)基于Windows操作系统且具有可视化界面的Earthview/InSAR、E砒)AS瓜ada饥nSAR和ENVI/SARScape。瑞士的GAMMA软件是一款综合性的商业InSAR数据处理软件,它能够实现从合成孔径雷达复数影像生成DEM、地表形变图以及土地利用分类图等数字产品的整个过程【521。GAMMA主要有两种相位解缠算法,分别是枝切区域增长算法以及最小成本流(MCF)技术与不规则三角网(唧)相结合的解缠算法。.Doris软件是由荷兰Delft大学Kampes等人使用面向对象的C++语言编写的,主 西南交通大学硕士研究生学位论文第6页要用来研究三维地形和地表形变【541。它利用第三方软件SNAPHU采用统计费用网络流算法进行相位解缠。由美国喷气推进实验室与加利福尼亚理工学院联合开发的开源软包Roi-pac(RepeatOrbitInterferometryPackage)也可进行InSAR数据处理,使用的相位解缠算法是经典Goldstein“枝切"法。EarthView/InSAR是加拿大渥太华AtlantisScientific公司开发的商业化InSAR数据处理软件,利用它可获得大区域甚至是全球范围内高精度的DEM和形变图,该软件使用的是迭代掩膜算法(IDM)来解缠相位【55】。美国ERDAS公司开发的遥感图像处理系统软件ERDAS脚舭SAR模块采用了最小费用流的解缠算法,该软件主要用来提取DEM。SARscape是由Sarmap公司研发的架构于专业遥感图像处理软件ENVI之上的雷达图像处理软件,提供了图形化操作界面,具有专业雷达图像处理与分析功能。它所包含的InSARJD.InSAR模块是通过区域增长算法或最小费用流算法进行相位解缠的。此外,2012年刘、黄分别基于路径跟踪和最小二乘相位解缠算法实现了相应相位解缠系统的开发[50,56]。伴随着合成孔径雷达干涉测量技术的发展,有关相位解缠算法和相位解缠软件也都得到了较为全面的发展,但是仍存在着不足之处。解缠算法种类多但普适性差,并且没有一个精度评估的标准,而InSAR数据处理软件所提供的的解缠算法均较单一,交互式操作也较少,不利于发现解缠过程中可能出现的错误并进行及时修正。1.3论文研究内容前已述及,相位解缠是InSAR和D.InSAR测量数据处理过程中的关键步骤与技术难点。相位解缠的准确性直接影响InSAR技术提取的地面高程以及D.InSAR技术测量的地表形变量的精度。针对上述对InSAR相位解缠研究及应用现状的分析,本文拟基于路径跟踪与最小范数相位解缠算法,开发一个专门针对相位解缠的综合应用系统,并为其添加解缠评估、DEM生成、图像处理、手动修改“枝切”以及框选解缠模块,将有助于提高解缠相位的精度,并可实现对相位解缠算法的综合评估,这将具有非常重要的应用价值与现实意义。具体来说,论文研究的内容主要包括以下几个方面: 西南交通大学硕士研究生学位论文第7页(1)在分析相位解缠原理的基础上重点阐述了四种路径跟踪相位解缠算法和四种最小范数相位解缠算法的基本原理与实现方法,并分析总结了不同解缠算法各自的特点、适用情形及优缺点。(2)在Windows操作系统下,基于面向对象的C群语言和GIS的嵌入式组件库ArcEngine,对InSAR相位解缠算法进行开发,并形成软件系统。(3)为提高相位解缠的成功率与可靠度,在相位解缠过程中引入人机交互式操作和评估策略,主要包括为系统添加解缠评估、DEM生成、手动修改解缠路径(即“枝切”)以及框选解缠模块。此外,为提高解缠效果,系统还添加了图像处理模块。(4)以各种模拟相位数据和真实干涉相位数据为数据源,对InSAR相位解缠系统中各种算法进行测试,通过对比分析验证算法的精度和可靠性。并且,系统将提供相位到高程的转换及地理编码功能,基于系统中各算法与SARscape中最小费用流算法的解缠相位提取出相应的DEM,通过二者的对比进一步验证解缠结果的精度,进而对解缠算法及所开发的系统的有效性进行评价。1.4论文组织结构论文的组织结构安排如下:第一章从InSAR和D.InSAR干涉测量理论入手探讨了相位解缠的研究意义,简要介绍了相位解缠算法和相位解缠软件的发展现状,并在此基础上陈述了论文的研究内容。第二章首先分别对基于路径跟踪的相位解缠算法和基于最小范数的相位解缠算法的基本原理进行阐述,并分析了与路径跟踪相位解缠紧密相关的留数点的概念,及留数点探测的方法。然后介绍了四种引导相位解缠的质量图和三种经典噪声滤波方法,最后还介绍了几种解缠评价的标准以及由解缠相位生成DEM的两个过程。第三章详细介绍了四种基于路径跟踪和四种基于最小范数的相位解缠算法的基本处理流程与实现方法,并对核心算法的实现进行了设计。第四章作为整个研究的核心,详细介绍了基于C撑和ArcEninge开发的InSAR相位解缠系统,并展示了各个模块的功能及相应的界面。此外,还对系统的框选解缠与手动修改“枝切”的具体操作流程进行了详细说明。 西南交通大学硕士研究生学位论文第8页第五章借助InSAR相位解缠系统针对模拟相位数据进行相位解缠实验,并详细分析了实验结果。第六章基于系统各算法的解缠相位以及系统提供的DEM生成模块所提取的DEM,并以SARsacape的解缠相位及其DEM为参考基准,评定了系统提供的八种解缠算法的精度,总结了各解缠算法的优缺点。最后是对论文所做工作的总结,并指出目前工作中所存在的一些问题及后续的研究方向。 西南交通大学硕士研究生学位论文第9页2.1引言第2章相位解缠原理及关键因素分析通过干涉或差分干涉所获取的干涉图中的相位值并非真实的绝对相位值,绝对相位值总是以非线性的方式被缠绕进(-万,万】之间,形成相位主值,即缠绕相位。缠绕相位公式如下:y(f)=缈O)+2万后(f)(2·1)其中,七(f)是一个整数函数,使一万<∥≤万;{f,(f)、缈(f)分别是缠绕相位和绝对相位。因此,为了得到真实的绝对相位值必须在缠绕相位值的基础上加上或者减去2万的整数倍,这个过程称为相位解缠(Ph嬲eUnwrapping)t571。由于一维相位解缠的积分路径唯一,并且是逐个积分,若受到相位噪声的影响或是遇到地形起伏导致相位不一致或不连续,使其中一点的解缠出现误差,那么这个误差会沿着积分路径传播到整个图像中【58,591。因此,为了寻求更多的路径来绕过不连续点,一维相位解缠逐渐发展为二维相位解缠。二维相位解缠需兼顾两个方面:一致性与精确性。一致性是指解缠后的相位数据矩阵中任意两点之间的相位差与这两点之间的路径无关;精确性是指解缠后的相位数据能真实地恢复原始的相位信息唧删。目前,常用的I.nSAR相位解缠方法可分为两类,即基于路径跟踪的相位解缠算法和基于最小范数的相位解缠算法。2.2基于路径跟踪的相位解缠基本原理基于路径跟踪的相位解缠算法是通过对相邻象元的缠绕相位梯度进行积分来恢复相位的真实值。二维相位解缠的线性积分方程可表示为:伊(,.)=LV伊·dr+9(%)(2-2)其中,9(,.)是点厂的绝对相位值,o,(ro)是起始点%的绝对相位值,C是区域D内连接%与,的任意路径,V9是9的相位梯度。根据微积分学: 西南交通大学硕士研究生学位论文第10页I=IF(r)·dr(2—3)式(2.3)的线性积分结果不仅依赖于路径C的起点和终点,还与路径C本身相关。如果函数F(,)在其定义域内处处可微,或旋度为0,等效为一个无旋场,那么该积分则与路径无关[63】。也即是说,在一条简单闭合路径积分,其结果应为零。JFCr)·dr=0(2-4),二维相位解缠一般采用上式来判断积分是否与路径无关。然而在干涉相位中,并不是任何积分都能够满足(2-4)的条件。有些像元上的缠绕相位数据由于噪声的影响或地形起伏比较剧烈,导致二维相位解缠并不孤立于积分路径。1987年Ghiglia、Mastin和Remero发现在二维曲面上沿着某一路径对缠绕相位差求积得到了不一致值(闭合路径积分的非零点被标记为不一致),解缠相位依赖于解缠路径,并且还发现不一致现象局限在某些孤立点或者孤立区域【641。后来,Goldstein、Zebker和Wemer在1988年将这些不一致定义为术语“留数点”,指出任何包含单个留数点或数量不相等的正负留数点的闭合路径积分都不为零∞51。不一致也被Huntley称为不连续[36,64]。由于留数点的存在,导致二维相位解缠依赖于积分路径,不同的相位积分路径导致不同的解缠结果,因而不能保证相位解缠的一致性与精确性。这也使得二维相位解缠成为一个相当难的问题,各种相位解缠算法其实都是在设法处理留数点,消除留数点对相位解缠的影响,从而获取精确的解缠相位。2.2.1留数定理二维相位解缠的留数定理:Ghiglia与Pritt将复变函数的留数定理应用N-维相位解缠中,则二维相位解缠的留数定理可表示成‘59】:【V妒(,)=2万×(闭合路径所包围的留数点电荷之和)(2.5)其中,留数点“电荷"表示一个整数,也就是说,围绕某一个留数点的闭合路径积分等于2万的整数倍。除病态相位结构外,留数点的电荷均等于±1。当积分区域内的留数点电荷平衡时,围绕该区域任意简单闭合路径的线性积分都为零,【V时)=o。因此,当且仅当所有积分路径都不包含未平衡留数点时,可得到符合一致性的解缠相位。利用“枝切"线将极性相反的留数点连接起来,使留数点电荷平衡,这样就能够 切"线(没有辅助信息或者假设条件),也不可能得知相位解缠具体在哪里出现问题(积分路径包围了未平衡留数点时解缠相位会出现不一致)。所以有必要研究留数点探测,2.2.2留数点探测上吡上0jj2A.3△ljr△3T卜:言△:···。o·1。。_Iji:。呻卜·o·2。_-△--:"◆0·4···假设有一个缠绕相位函数矩阵杪(朋,疗),图2-1显示的是Nj}gW-d,N#。图中所示的是经归一化处理后的缠绕相位值,在(-0.5,0.512.N,乘以2万才是以弧度为单位的Al=形缈(聊,刀+1)一缈(所,刀))A22矿{少(研+l,船+1)一y(m,刀+1)>(2.6)A3=形{y(m+1,刀)一缈(,行+1,,l+1)>、7A4=∥{y(肌,n)-Ic(m+l,刀))然后,对最小闭合路径(实线箭头所示)的缠绕相位梯度值求和计算留数点电荷。 西南交通大学硕士研究生学位论文第12页4g=∑△,=△l+△2+△3+△4=-0.2-0.1+0.4-0.1=0(2—7)f;l此例中,留数点电荷为0,表示相位是连续的,不存在不一致现象,该闭合路径的左上角像素不是留数点。再计算虚线箭头所示路径的留数点电荷,可以得到:4q=∑A:=△i+△:+△:+△:=-0.4-0.2-0.3-0.1=一1.0(2-8)1=1在此找到了一个非零积分路径,将该积分路径的左上角像素标记为留数点,并且标记其电荷极性为负。如果在整个缠绕相位矩阵y伽,刀)上进行解缠运算时,这种由简单闭合路径积分所造成的误差会依次传播到整个矩阵中去,从而造成整个图像的误差。因此,若能正确检测出所有留数点,并用“枝切”线将正负极性的留数点连接起来,在保证不穿过“枝切"线的情况下,可以沿任意路径进行积分,并且得到的结果相同。“枝切’’线是以正(负)极性留数点开始,以负(正)极性留数点终止(除遇到相位数据边界外)。根据这一原则可先建立相位数据矩阵的“枝切”线,并要求任一“枝切”线的正负极性留数点的距离最近,即满足“枝切”线最短准则。进一步分析得知:任何积分路径如果包含不同个数的正负留数点,即该闭合路径积分区域内的留数点电荷未达到平衡,则沿该路径解缠的结果不能满足缠绕相位梯度连续的要求;当正极性和负极性留数点个数一致时,则可以得到连续正确的解缠结果。2.3基于最小范数的相位解缠基本原理最小范数解缠算法完全不同于路径跟踪解缠算法,它着眼于整体,采用最优化思想,寻求最小二乘意义下的最优解缠结果。该算法不用识别留数点,它的基本思想是使缠绕相位梯度与解缠相位梯度在x和Y方向上的差达到极小。假设给定二维缠绕相位函数帆。,(扛0,1,...,M-1;j=0,1,...,N-1),其对应的解缠相位函数为谚。,O=o,l,...,M-I;j=o,1,...,N-1),则最小r范数解满足:M一2Ⅳ-1M-1Ⅳ-2Min{J}=Min{∑I谚+。√一谚,,_AuXIp+∑∑协,p.一谚,,一△0l,)(2-9)i=Oj=OI=0j=O 西南交通大学硕士研究生学位论文第13页通过对上式中谚,.,求偏微分并经整理可得到最小范数问题的广义解:其中:(谚+l,,一识,_,一△i,)uq,_,)+(谚,,+l一谚,,一△乙)矿G,jf)一(谚√一谚一l√一△:l√)u(f一1,jf)一(谚√一谚√一l一△01)y(f,_,一1)=0(2-10)【,(L_,)={10谚+l√一谚√一△。r一207t=he。r’1w”is”e’M一2;/=。,1’⋯,Ⅳ一1(2.11)l、-一7y(毛/):{I力√+-一谚√_Aby,lp一2‘=o,1””,M-1;jf=o,1,⋯,Ⅳ一2(2.12)10otherwise当P=2时,可由方程(2—10)得到泊松方程的离散形式:其中:(谚+l,J一2办,J+谚一l,/)+(谚,J+1—2谚。.,+谚,产1)=房,.,(2—13)岛,』=(△i,一Al。_ta)+(△乙一△乙一1)(2一14)由此可见,当P=2时,可根据最小二乘准则将二维相位解缠问题转化为求解纽曼边界条件下的离散泊松方程。无权重最小二乘相位解缠算法是通过使缠绕相位梯度与真实相位梯度的差的平方和最小来实现相位解缠,将局部噪声对整个解缠平面的影响降到最小。M-2N-1M-1N-2s2=∑∑(谚+l√一谚,』一△0)2+∑∑(谚,_,+1一谚√-Ay“)2(2—15)i=Oj=oi=0j-o因此,任何可以有效求解泊松方程的算法都可以用于相位解缠。然而无权重最小二乘相位解缠算法虽然运算速度快,但由于解缠时没有绕过而是穿过相位数据中不连续区域,因此造成局部误差的传播,从而引起全局误差。这个缺陷可以通过引入权重来弥补。一般说来,权重数据Ⅵ。,的取值范围是【0,1】,权重由相位质量图或者其他先验知识来确定,则最小二乘问题可转化为加权最小二乘问题:^f一2^r—l^,-1N-2s2=∑∑U(i,_,)(谚+l,J一谚。/一△i』)2+∑∑矿o,歹)(谚,p。一谚。,一△乙)2(2·16)i--Oj=Oi=Oj=0 西南交通大学硕士研究生学位论文第14页其中,梯度权重u(i,J)、g(i,J)定义为:Uw(i棚,y);--蔬彩W2陋∽y(f,歹)=n曲(幢,+l,u)这个问题的最小二乘解谚。,定义如下:u(‘,/)(谚:-。.,一谚。?)一uO一1,,)(谚?一谚q√)(2.18)+矿(f,j『)(谚,.,+l一识√)一y(f,_,一1)(谚,/一绣,/-1)=q√、’式(2-18)中,ct。,为加权相位拉普拉斯操作数:q,/=u(f,j『)A,xj—U(i一1,歹)A二l,,+y(f,歹)△乙一矿(f,J一1)△i,一l(2—19)对式(2—19)进行整理得:破,:堕丛型竺兰堕趔趔二!监.墨:2(2-20)“”u(i,j)+U(i-1,jf)+V(i,j『)+V(i,J一1)、7对式(2-20)进行多次迭代可求得力。,。加权最小二乘问题还可用矩阵形式表示成:刨=C(2-21)2.4关键因素分析2.4.1质量图在二维相位解缠中,质量图是一个很重要的概念,很多相位解缠算法都利用质量图来指导解缠,从而提高解缠的准确性。如质量图路径引导算法以及加权r范数算法等。质量图是用来评价相位数据质量的数组,定义了相位图中每个像素对应相位质量的好坏。常用的质量图主要有四种:相干系数图、伪相干系数图、相位导数变化图以及最大相位梯度图。除相干系数图外,其余三种质量图均可利用相位解缠系统从干涉相位中直接提取。相干系数是由两幅SAR影像定义,表明不同区域的相干性,是干涉质量最客观的评价,可作为质量图的首选。通过各种实验表明:在没有相干图的情况下,当地形起伏较大时,相位导数变化图可作为质量图第一选择,其次分别是最大相位梯度图、伪相干系数图;当地形较为平坦时,伪相干系数准确反映出了噪声的位置,质量图应首先选择 西南交通大学硕士研究生学位论文第15页伪相干系数图,然后依次是相位导数变化图、最大相位梯度图‘66,671。2.4.2掩膜在相位解缠过程中,有的区域由于失相关严重(例如,水面等),导致相位不连续;有的区域地势平坦不需要进行滤波。在这样的情况下,可以使用一个掩膜图将这些区域掩盖起来,让这些区域不参与解缠和滤波,以防止误差的传播。掩膜图是利用质量图来制作的,是将质量指标二值化为“O”和“1”的数组。“0’’表示低质量像素或零权值点,即被掩盖掉的不再参与解缠的像素;“1"表示高质量像素,即未被掩盖的像素。掩膜图的生成方法一般采用阈值法,通过对质量图给定一定阈值进行判断从而得到一个掩膜:像素质量值小于阈值的被标记为“0";像素质量值大于阂值的就被标记为“1”。基于此,InSAR相位解缠系统在提取相位质量图时也伴随生成了相应的掩膜。2.4.3噪声滤波干涉图中如果存在噪声,就有可能造成相位数据的不连续性和不一致性,使得相位解缠过程中产生局部误差,并在解缠过程中传播,对解缠精度和效率产生很大影响【吲。而噪声的存在常常是不可避免的,因此在生成相位干涉图后,一个重要的工作就是尽可能滤除噪声,提高信噪比,减少留数点,以提高相位解缠的精确性与一致性嗍。干涉图噪声的滤除要在提取相位值之前的复数中进行,不能直接处理相位值【701。如果直接处理相位值,干涉条纹边缘会被平滑,就不能保持较好的陡峭形状。因此,下面介绍的三种经典滤波算法都是在干涉相位复数数据上进行的。(1)均值滤波均值滤波是在干涉相位上选取一个kxk大小的窗口,计算该窗口内所有像素复数值的实部之和、虚部之和,再对其求反正切即为该窗口中心像素(m,刀)的相位值:m+gn+g册+七n+lt多册一=arctan(E∑sill%∥∑∑cos%√)(2.22)f-历j-ni--m』=”式中,sin5f6。.,表示点(f,.j『)复数值的虚部,cos%,.,表示其实部。(2)自适应滤波自适应滤波只选取了取样窗口内的某些像素参与平均计算,例如Sigma滤波。这 西南交通大学硕士研究生学位论文第16页种滤波只适用于InSAR数据,因为它需要利用相干系数来定义哪些像素参与平均计算。给定一个相位(myn)的相干系数7,计算以(脚,刀)为中心的kxk窗口内所有像素相干系数的标准偏差仃,然后选取该窗口内相干系数位于(y-20",y+20")之间的那些像素的相位值来进行平均计算。对Si印m滤波进行改进,利用其它质量图(例如:相位导数变化图、最大相位梯度图、伪相干系数图等)代替相干系数,这样的自适应滤波可适用于任何相位数据。相位解缠系统采用的即是这类滤波。(3)中值滤波当相位数据存在“椒盐’’噪声时,中值滤波比均值滤波更容易减弱噪声。由于复数没有一个自然的排序,因而没有一种自然的方式来定义中值滤波。七×七窗口内中心像素(m,刀)在中值滤波中相位值计算的近似方程为:y用JI=arctan(median{sinqzl.i},median{cosq/i,i))(2-23)与均值滤波相比,中值滤波能够较好的保护相位条纹,同时能够有效滤除颗粒噪声,但是没有利用信号的统计特性,所以结果并非最优【701。为了使相位解缠更加容易,解缠效率更高,相位解缠系统的图像处理模块也提供了自适应滤波来对解缠前的干涉图进行去噪处理。可根据需要滤波的干涉数据具体情况,选择合适的相位质量图来进行滤波。2.5解缠评价标准解缠算法优劣的评定应包括其精度与准确度、计算效率、内存使用等情况的考察。本文主要采用解缠误差图、不连续图和反缠绕相位图三种指标来评价解缠结果。2.5.1误差图对于模拟数据来说,可以直接利用解缠相位与参考绝对相位做差生成的误差图来评估解缠结果。但是对于实测数据无法获得绝对相位,这时可以将Goldstein“枝切"算法的解缠结果作为基准,将其与其他算法的解缠结果进行求差生成误差图,以此来分析不同算法之间的差异。这是由于在干涉相位数据不连续区域较少的情况下,Goldstein“枝切"算法的解缠结果普遍被认为是最好的【32删。因此,本文采用了两种相位数据求差的误差图作为评价标准,以此来分析比较不同解缠算法的优劣。 西南交通大学硕士研究生学位论文第17页2.5.2不连续图Goldsten“枝切”算法中生成的“枝切’’线可以用于确定相位解缠中出错的部位,而不连续图就在解缠结果评价中起着“枝切"分布图的作用,它指出了相位解缠路径未经过的区域,通过不连续图可以定性的看出解缠结果的好坏。并且由于“枝切"分布图只在Golds'ten“枝切’’算法和“掩膜枝切"算法中才能得到,而不连续图可利用任何算法的解缠相位计算生成,因此本文采用了不连续图作为评价解缠算法的标准。2.5.3反缠绕相位图反缠绕相位图可从视觉上定性地说明解缠结果的准确性,同时也能说明最小二乘算法对无噪声相位数据的缠绕方式。反缠绕相位图一般用于最小二乘算法,而对路径跟踪算法则无多大意义。因为最小二乘解与初始数据不一定全等,如果两者存在差异,那么反缠绕相位图会将这些不一致表现出来。但是,路径跟踪解的反缠绕结果与初始值是一致的,所以一般不对路径跟踪算法的解缠结果使用反缠绕相位图来评价。2.6InSAR解缠相位的应用—-DEM的生成常规InSAR相位解缠结果可用于数字高程模型(DEM)的生成,并且,相位解缠的结果直接影响到生成的DEM精度,那么可以通过生成的DEM的精度来评价相位解缠的准确性。因此,相位解缠系统专门添加了DEM生成模块,一方面用来评估相位解缠的精度,另一方面也方便用户直接得到最终的InSAR数据处理产品——数字高程模型(DEM)。DEM的生成主要分为两个子模块完成:生成高程图和地理编码。2.6.1生成高程图根据式(2.24)相位与地面高程的关系,可实现解缠相位值到地面高程值的转换,但是此时得到的DEM仍是斜距坐标下的,还需要将高程数据从斜距坐标转换到地图坐标系统,在此过程中要进行几何校正、数据重采样,才能得到DEM。办=等警陋24,其中,矽为绝对相位;允为雷达波长:秒为雷达侧视角;毋为基线在斜距墨上的垂直投影分量。 西南交通大学硕士研究生学位论文第18页2.6.2地理编码地理编码即是对高程数据应用几何变换,将它从雷达几何结构转换为地图结构,再对其进行重采样得到与地形图相符的DEM。主要包括两个步骤:(1)利用构像方程解决定位问题给定图像点的行列号坐标(,,,)和高程h,根据构像模型,求解该点对应的地面坐标。本文采用了距离一多普勒(Ranger-Doppler,RD)模型来建立影像坐标与地面坐标间的数学关系,即影像与地面间的坐标关系,然后利用牛顿迭代法解算RD模型,最后将迭代计算得到的地面点的地心坐标转换为大地经纬度坐标。具体的距离.多普勒构像模型【7l】如下:1)地球椭球方程器+丽Z2一协25,(口+办)2(6+办)2、7其中,a为地球椭球长半轴;b为地球椭球短半轴;h为地面点高程。2)多普勒方程血一争错陋26,其中,厶为雷达波束的多普勒位移;五为雷达波长;足、圪分别为卫星的位置矢量和速度矢量;鬲、露分别为地物点的位置矢量和速度矢量。3)距离方程R=I可一面J-√(K一置)2+(鬈一z)2+(五一z|)2(2.27)(2)高程重采样将原始图像中某像素的高程值转化为DEM中该像素对应位置处的新高程值的过程,一般采用内插算法,如双线性内插、三次卷积内插等。但是由于内插所需的点较难确定,所以相位解缠系统选用GIS中的不规则三角网(TriangulatedIrregularNetwork,TIN)来生成DEM。TIN采用离散数据点生成连续不重叠的不规则三角网格来表示地形表面,在地形平坦区域三角形较少,而在地形复杂的区域三角形较多‘731。因此,TIN 西南交通大学硕士研究生学位论文第19页能较好的顾及地形地貌特征,逼真表示复杂地形的高低起伏变化,并且能克服地形平坦区域的数据冗余【74】。2.7小结本章分别阐述了基于路径跟踪相位解缠算法与基于最小范数相位解缠算法的基本原理,并分析了与路径跟踪相位解缠紧密相关的留数点的概念及留数点探测的方法。然后介绍了引导相位解缠的四种质量图,利用这些质量图还可制作掩膜,将干涉图中相位不连续的区域掩盖掉,让这些区域不参与解缠。为了使相位解缠更加容易,解缠效率更高,本章还介绍了三种经典的噪声滤波方法。针对目前相位解缠仍没有一个精度评估标准的情况,本章提出了三种解缠评价的标准,分别采用解缠误差图、不连续图和反缠绕相位图三种指标来评价解缠结果。此外,为了进一步分析各种解缠算法的精度,本文还采用了根据各算法的解缠相位提取出的DEM的精度来评价算法的准确性。为此,解缠系统还专门添加了一个DEM生成模块。因此,在本章最后详细阐述了解缠相位生成DEM的两个过程。 西南交通大学硕士研究生学位论文第20页3.1引言第3章相位解缠算法目前,国内外学者已研究出了众多的相位解缠算法,按解缠方式的不同大致可分为两类:基于路径跟踪的相位解缠算法和基于最小范数的相位解缠算法。本章将对这两类方法的代表性算法的处理流程和实现方法进行详细的阐述。3.2基于路径跟踪的相位解缠算法基于路径跟踪相位解缠算法的基本策略是将可能的误差传播限制在噪声区域,通过识别并连接留数点以选择合适的积分路径,隔离噪声区,防止相位噪声的全程传播。主要的算法有:Goldstein“枝切"算法、质量图路径引导算法、“掩膜枝切"算法、Flyrm最小不连续算法等。3.2.1Goldstein“枝切’’算法Goldstein“枝切"(Goldstein'sBranchCut)算法的基本思想是:识别正负留数点,并以最近邻的方式连接极性相反的留数点对或连接包含多对极性相反的留数点对的片区,实现留数点电荷的平衡,生成最优“枝切"线,确定积分路径,防止误差沿积分路径传播。此外,也可通过连接留数点与相位数据边界的方式,实现留数点电荷极性的抵消。Goldstein“枝切”算法基本步骤:Stepl.识别留数点;Step2.生成“枝切”线;Step3.围绕“枝切”线进行路径积分。Goldstein“枝切"算法实现方法的相应伪代码如图3.1至3-4所示。‘'UpdatetheAdjoinList'’子程序利用某个像素作为输入参数,并判断相邻的四个像素是否应当被放入AdjoinList中以待解缠,相应的伪代码如图3.1所示. 西南交通大学硕士研究生学位论文第21页图3-IUpdatetheAdjoinList的伪代码图3.2围绕“枝切”线进行路径积分的伪代码 西南交通大学硕士研究生学位论文第22页for(每一个未被平衡的留数点){将该留数点标记为ACTn吧(该留数点被标记为POS—RES)7charge=1:charge=一1for(N=3至最大的窗口大小){for(每个被标记为ACTⅣE的留数点){以当前活动留数点为中心安置N×N大小的窗口for(N×N窗口内的每个像素:窗口像素){if(该窗口像素被标记为BORDER){il:charge=0在当前活动留数点与该边界像素之间安置“枝切”)elseif(该窗口像素为一个非活动留教点){if(该窗口像素未被平衡){将该窗口像素的极性累加至charge:(该窗口像素点带POS—RES#{g)?charge+=1:charge+=一1将该窗口像素标记为ⅥSITEDl将该窗口像素标记为ACTⅣE在当前活动留数点与该窗口像素之间安置“枝切”.’if(charge=0)跳至下面标记·处)*if(charge不等于o)在未被平衡的留数点与边界之间安置”枝切”为所有带ACTIVE标记的像素加上VISITED标记并取消ACTIVE标记}图3-3“枝切”线生成的伪代码 西南交通大学硕士研究生学位论文第23页for(干涉相位数据中的每个像素(朋,刀)){对最小闭合路径内的缠绕相位梯度值求和(形表示缠绕操作):W{g/(m,刀+1)一y(埘,刀)}+形{y(所+1,刀+1)一∥(所,刀+1))+形{y(,”+1,刀)一g(m+l,刀+1))+W{g/(m,万)一g(m+l,刀)’if(和为2万)将该像素(聊,丹)标记为POS—RESif(和为-2zt)将该像素(埘,刀)标记为NEG—RES,3.2.2质量图路径引导算法图3_4留数点识别的伪代码路径跟踪算法除了通过检测留数点的分布来确定积分路径外,还可利用一些辅助信息来分析相位数据的质量,更好的设置积分路径,相位质量图就是一种很好的辅助信息。质量图路径引导算法(Quality-GuidedPathFollowing)就是利用质量图来指导积分路径的生成,使得积分路径总是沿着高质量像素开始,避开低质量像素,最大限度地阻止误差传播。质量图路径引导算法基本实现过程如图3.5至图3.7所示。图3-5UpdatetheAdjoinList的伪代码 西南交通大学硕士研究生学位论文第24页图3-6质量图路径引导算法实现的伪代码图3-7InsertPixelinAdjoinList的伪代码3.2.3“掩膜枝切”算法“掩膜枝切”(Maskcut)算法结合了Goldstein“枝切”算法和质量图路径引导算法各自的优势:识别留数点并用“枝切’’线将留数点连接起来,防止解缠路径包围不 西南交通大学硕士研究生学位论文第25页平衡的留数点;利用质量图来有效引导“枝切"线的设置。该算法既充分利用了外部辅助信息,同时也确保了解缠相位不因积分路径包围留数点而产生错误的2a整数倍的跳变。“掩膜枝切"算法基本步骤:Stepl.识别留数点;Step2.生成“掩膜枝切”;Step3.细化“掩膜枝切";Step4.围绕“掩膜枝切”进行路径积分。“掩膜枝切’’算法实现方法的相应伪代码如下:for(每个未被平衡的留数点){将该留数点标记为VISITED以表示已被平衡(该留数点的带POS—I也S标记)7charge=1:charge=-1将该留数点标记为BRANCHCUT(即为maskcut)更新adjoinlist(见图3-10的伪代码)while(a两oinlist不为空&&charge≠0){从喇oinlist中获取最低质量像素将该像素标记为BRANCHCUTif(该像素为一个为平衡的留数点){(该像素带POS—I也S标记)7charge+=1:charge+=-1为该像素添加标记VISITED以表示已被平衡'if(该像素是一个边界像素)charge=0更新删oinlist(见图3-10的伪代码))图3-8生成“掩膜枝切”(maskcut)的伪代码“掩膜枝切骨是利用区域增长法生成的,素。为了生成优化程度较高的“掩膜枝切",像素的“枝切一标记【501。比较厚,其中包含了许多非“枝切"像“ThintheMaskCuts”程序将取消这些 西南交通大学硕士研究生学位论文第26页图3-9细化“掩膜枝切”(ThintheMaskCuts)的伪代码图3一10UpdatetheAdjoinList的伪代码3.2.4Flynn最小不连续算法缠绕相位数据(干涉图)表现为一系列的“条纹线”,其中最亮像素(假彩色中的红色)与最暗像素(假彩色中的蓝色)之间的“条纹线"标志着一万到万的转换,这些“条纹线’’将图像分割成若干区域,形成一系列的“不连续"。这样,相位解缠的过程也可以看作是为每一个分割区域分配2万整数倍,使相位的“不连续”性达到最小。1997年Flynn提出了Flynn最小不连续算法(Flyrm’sminimumDiscontinuity),采用树节点的方法识别“不连续线”而不是“条纹线",探测“不连续性”连成的“圈"(100p),通过迭代计算为“圈"所包围的相位添加适当的27r整数倍,达到最小化“不连续"的目的【39】。 西南交通大学硕士研究生学位论文第27页Flyrm算法基本步骤:Stepl.根据输入的缠绕相位数据计算跳变数;Step2.重复扫描节点,添加“边”并移除“圈",直到没有新的“边’’被添加也没有“圈"被移除;Step3.根据跳变数计算解缠相位。Flynn算法具体实现方法相应的伪代码如图3-11至3-14所示。图3-11RemoveLoop(base,last)的伪代码图3·12Cllang嘶删似丘口m,t0)的伪代码 西南交通大学硕士研究生学位论文第28页设置所有”节点值”=0设置所有相邻的”节点对(from,t0)”的边标记edge(6∞m,t0)=falserepeat{设置“边”发现标记为edge—found=faslefor(所有相邻的“节点对(from,t0)”){if(该”节点对”的边标记edge(右的m,to)=false){b'V(from,to)=node(from)0node(to)的“边(edge)”值的累加值纵节点”到“到节点”的“节点”值变化量:dval=”从节点”程t+SV(from,t0)·”到节点”值if(dval>0){edge—found=true设置”圈”发现标记100p—found=ChangeValue(to,dval,from)(调整整个“节点”数组的值)if(100p—found=true){RemoveLoop(to,from)(除已发现的围成“圈”的所有边界,并更改相应的跳变数)ChangeValue(t0,-dval,Null))else{for(除”从节点”以外的“到节点”的其余三个相邻像素(用nbr表示)){边标记e趣e(nbr,to)=false(即取消点nbr到,,从节点”的边标记,表明无边存在)edge(from,to)=true)1)’)until(不再发现新的边界:edge—found=false)图3-13Flynn最小不连续算法步骤二的伪代码 西南交通大学硕士研究生学位论文第29页图3.14ChangeValue(base,dval,100p)的伪代码3.3基于最小范数的相位解缠算法最小范数解缠算法着眼于整体,采用最优化思想,寻求最小二乘意义下的最优解缠结果嘲,主要包括无权重最小二乘相位解缠算法、加权最小二乘相位解缠算法以及最小r范数相位解缠算法等。3.3.1无权重最,b--乘算法由于SAR影像维数较大,对于常用的迭代算法,如Jaeobi迭代、Gauss.Seidel迭代、SOR迭代等,其求解效率较低,不适合实际数据的计算。下面主要介绍两种高效率的相位解缠算法:基于离散余弦变换∞CT)的算法和无权重多重网格算法。(1)离散余弦变换算法对二维离散函数谚,,进行余弦变换(DiscreteCosmeTransformation,DCT),得到如下结果: 西南交通大学硕士研究生学位论文第30页C埘一=_a篓t1-1势Jv-A严s[刍酢㈣cos[蠢酢删善丢4谚√c。s[奇酢川)】c。s[紊酢川)】0≤m≤M一1,0≤珂≤N一1(3.1)0otherwiseq。。进行二维DCT逆变换后,其结果谚,J为:谚,J=面1M刍-IN刍-1w(朋,刀)%cos【刍肌(2⋯)】c。s爵犯川)】0≤i
此文档下载收益归作者所有