毕业设计（论文）MODIS 1B影像的BowTie效应快速去除技术

《毕业设计（论文）MODIS 1B影像的BowTie效应快速去除技术》由会员分享，可在线阅读，更多相关《毕业设计（论文）MODIS 1B影像的BowTie效应快速去除技术（30页珍藏版）》请在装配图网上搜索。

1、MODIS 1B影像的Bow-Tie效应快速去除技术Rapid Removal of Bow-Tie Effects from MODIS 1B Data: An Algorithmic Development目录摘要-3Abstract-4一、 引言-5二、 算法设计-6（一） Bow-Tie现象成因分析-6（二） MODIS 1B的纬度科学数据集-7（三） Bow-Tie去除算法流程-8三、 算法实现-12（一） 纬度数据的蝴蝶结现象去除代码示例-12（二） 图像数据的蝴蝶结现象去除代码示例-14四、 结果与分析-15五、 结论与讨论-17参考文献-19致谢-19附录1程序运行界面-20附

2、录2程序源代码-21MODIS 1B影像的Bow-Tie效应快速去除技术摘 要：中分辨率成像光谱仪（MODIS）是一款广泛应用于进行大、中尺度大气、海洋和陆地综合监测的新一代被动式摆动扫描遥感传感器，其55o的宽视场角和2330公里的幅宽使其影像边缘区域存在明显的重叠状几何畸变，称为Bow-Tie效应或蝴蝶结效应。目前该类几何畸变通常在进行投影转换时同步去除，但投影转换处理耗时较长，不利于影像的直接快速浏览。而现有Bow-tie效应去除算法除了ENVI外多直接基于影像上重叠区的检测，未充分利用MODIS自带的经纬度信息，算法效率不高。本文从MODIS影像的成像机制及其1B文件的存储结构入手，根

3、据其自带的星下点1公里间隔纬度数据计算扫描方向上各列的重叠程度，并以此为依据进行Bow-tie效应去除。此外，鉴于部分纬度数据可能由于仪器等原因出现异常值，在本文的算法中也对这些异常进行了纠正。本算法在Visual C+ 6.0环境下开发实现，数据处理速度与纠正后影像质量均达到了设计要求。关键词：MODIS; Bow-Tie效应去除; 算法设计; 畸变去除Rapid Removal of Bow-Tie Effects from MODIS 1B Data: An Algorithmic DevelopmentAbstract: Moderate Resolution Imaging Spec

4、troradiometer (MODIS) is a passive swing spectral scanner of new generation for collecting remote sensing data to be used in building a global dynamic model of the earth environment. Because of its 55o wide observation angles, obvious overlapping, known as the Bow-Tie effect, occurs in the edge area

5、 of MODIS scanned images, which are provided to end users as the 1B data. The existing algorithms for the removal of the Bow-Tie effect usually integrate this process into the image projection, which takes a considerable amount of time before seeing the results. In addition most of these algorithms

6、are based on direct detection of the overlap between adjacent swaps and do not utilize the relevant geolocation information contained in the HDF format. The study presented here attempts to develop a computer algorithm by focusing on the imaging mechanism of MODIS and the data structure of the 1B fi

7、le. In this algorithm, the Bow-Tie effect is removed according to the overlap of its 1KM latitude coordinate values. Furthermore, this algorithm also provides functions for missing data examination and correction. This algorithm was implemented using Visual C+6.0. Both the effectiveness and efficien

8、cy of the algorithm for removing the Bow-Tie effect were evaluated. The experiments demonstrated a satisfactory image quality after the removal processing and a significant improvement of execution time over some major algorithms for the same purpose.Keywords: MODIS; Bow-Tie effect; algorithm design

9、; distortion removal 一、 引言中分辨率成像光谱仪（MODIS）是由NASA联合美国和国际多个科研机构历时8年、耗资13亿美元研制的新一代星载遥感传感器，其目标主要是建立全球大气、海洋和陆地的动力学模型，并能在地球发生变化之前做出准确的预测。MODIS搭载在地球观测系统(EOS)的Terra和Aqua两颗卫星上 刘玉洁，杨忠东. MODIS遥感信息处理原理与算法M. 北京：科学出版社，2000.。EOS计划的目标是由一系列低轨道大型卫星平台装载先进仪器，以获得大量的先进卫星遥感资料，进而通过反演获得能反映地球这个复杂系统变异的多方面综合信息。MODIS共有36个观测通道，波

10、谱范围从0.4m到14m，覆盖范围包括可见光、近红外、短波红外和远红外谱区。MODIS的星下点分辨率可分为250m、500m和1000m三种。目前MODIS共有44种标准产品，按研究对象可分为陆地、海洋和大气三类；按接收和处理过程又可分为1A、1B、L2、L2G、L3和L4等多个级别 Modis Level 1A Earth Location: Algorithm Theoretical Basis DocumentE. 1997, http:/modis.gsfc.nasa.gov。其中经辐射定标和地理坐标定位的MODIS 1B数据是本文的研究数据源。华东师范大学地理信息科学教育部重点实验室

11、在自然地理馆楼顶安装了卫星数字广播接收系统(DVB-S)，可准实时接收国家卫星气象中心汇集处理后发布的MODIS 1B数据。本文所采用的示例数据为2006年10月9日国家卫星气象中心北京站接收的Aqua/MODIS数据。目前在国内外针对蝴蝶结现象去除的专门研究开展不多，绝大多数将蝴蝶结现象的去除与MODIS数据的几何纠正同时进行。这样虽然可直接利用现有几何纠正算法，如多项式纠正模型等，但是其执行效率却非常低。国外针对MODIS 1B数据的处理软件MRT(Free MODIS Reprojection Tool)以及通用遥感影像处理软件ENVI均可以进行这样的处理（注：ENVI 4.2版本后中的

12、相关算法已开始基于地理坐标进行处理，但仍未将蝴蝶结效应去除与图像及和纠正独立出来），但根据处理的数据分辨率不同，其耗费的时间往往需要数分钟甚至是几十分钟。对于单独的蝴蝶结现象去除的研究相对较少。郭广猛(2003)提出了利用IDL语言开辟一个33的滑动窗口，来对整幅影像进行扫描判断改矩阵的相关系数，并以此对重叠的图像进行纠正；程亮等(2005)利用MODIS 1B数据所对应的MOD03数据即地理修订文件来对图像进行扫描，并利用分形插值方法进行蝴蝶结现象的去除。其他的相关研究基本也是围绕着两种方法展开。这两种方法或者仅从图像学的角度进行蝴蝶结现象的去除或者利用MODIS的星历文件进行去除，算法设计

13、较为复杂，执行效率也相对较低。因此在MODIS数据接收终端等应用环境下，如何快速地对蝴蝶结现象进行去除，增强影像的可读性便成为本文的主要研究目的。二、 算法设计（一） Bow-Tie现象成因分析MODIS探测器是一种被动式摆动扫描探测器，其横向扫描角为55。每条扫描线由1 354个1 km的像素或2 708个500 m像素或5 416个250 m像素组成。由于地球曲率的影响,扫描线的实际跨度大约为2 330km。每完成一次扫描,MODIS探测器沿轨道前进10 km,这10 km的航向区域就构成了一个扫描条带。每个扫描条带的宽度分别为10个像素(1 km分辨率)、20个像素(500 m分辨率)和

14、40个像素(250 m分辨率) Running SW, etc. Terrestrial Remote Sensing Science and Algorithms Planned for EOS/MODIS J. INT. J. Remote Sensing, 1994 Vol. 15, No. 17, p3587-3620.。由于MODIS探测器对地球观测的视野几何特性、地球表面的曲率、地形起伏和MODIS探测器运动中的抖动等因素的共同影响, MODIS 1B数据存在多种几何畸变, 其中MODIS 1B数据相邻扫描带之间存在的错位和重叠现象（即蝴蝶结效应）是其重要的系统畸变之一。“蝴蝶结”

15、现象随着观测角度的增大而趋严重,计算表明,在视场角为24时,相邻条带像元间的重叠度为10%,而在条带的两端视场角为55时,重叠度可达50% Robert E Wolf, David P Roy, Eric Vermote. MODIS Land Data Storage, Gridding, and Compositing Methodology: Level 2 Grid J. IEEE Transactions on Geoscience and Remote Sensing, 1998, 36(4):p1324-1338.。其基本的形成原理如图1所示。图中显示有3个扫描带(横向),每一条

16、扫描行有10个像元,每个像元在星下点处分辨率为1 km,所以每个扫描行宽度为10km。由于每扫描一次地球会自转一定距离,图上显示3个扫描行之间在东西方向上有一个小的错位。由于地球曲率的影响和扫描角度的增大,越向边缘像元尺寸越大,在边缘处一个扫描带的宽度已增大到20 km,相邻的扫描带之间已有10 km的重叠 李登科. 消除MODIS图像重叠现象的方法研究J. 山西气象, 2005, (3), p1-3。图1MODIS 1B 1KM数据左半部分蝴蝶结效应示意图Figure 1. The Bow-Tie effect on left part of the MODIS 1B 1KM data通过上

17、述分析可见，正是由于MODIS的这种成像机制，决定了其影像不可避免地带有蝴蝶结效应，且该效应越往两侧越严重，这种系统性几何畸变大大影响了影像的可读性，尤其是影像边缘区的重叠现象几乎完全掩盖了地表的真实情况，对MODIS 1B数据蝴蝶结效应的去除是对MODIS数据进行预处理的必备操作。（二） MODIS 1B的科学数据集MODIS 1B数据主要以HDF文件的格式进行存储和分发。HDF（Hierarchical data format）格式是由美国伊利诺大学的NCSA（National Center for Supercomputing Applications）组织开发的一种超文本文件格式，是一

18、种用于存储和分发科学数据的自描述性、多对象的层次数据格式 景毅刚 等，HDF-EOS数据格式及其数据导入J, 气象科技，2005-6 333，p278-282。通常其记录的对象包含多维排列,用于存储科学数据，即称为科学数据集（SDS）。1B级使用SDS存储校准过的科学数据,它们的不确定性,部分质量保证数据。通过属性数据使SDSs能够自我描述。利用HDF所提供的相关库函数可以很方便将文件中存储的原始数据读入内存。由于MODIS 1B数据自带有经纬度信息，并以1KM分辨率存于1KM、500m、250m各个分辨率的数据文件中 Castleman K R. Digital Image Processi

19、ngM. Englewood Cliffs,NJ: Prentice-Hall, 1979.致谢：本论文在算法设计、程序编写时得到施润和老师、吴健平老师的指导与帮助，在此向两位老师的辛勤工作与无私帮助表示衷心的感谢！同时向审稿老师以及周围的同学给与的指导与帮助表示由衷的感谢！附录1程序运行界面附录2程序源代码/ FastBTRemoveView.cpp : implementation of the CFastBTRemoveView class/#include stdafx.h#include FastBTRemove.h#include iostream.h#include FastBT

20、RemoveDoc.h#include FastBTRemoveView.h#include afxcmn.h#include MODISHDF.h#include Afxwin.h#ifdef _DEBUG#define new DEBUG_NEW#undef THIS_FILEstatic char THIS_FILE = _FILE_;#endifextern MODISHDF theHDF;/ CFastBTRemoveViewIMPLEMENT_DYNCREATE(CFastBTRemoveView, CView)BEGIN_MESSAGE_MAP(CFastBTRemoveView

21、, CView)/AFX_MSG_MAP(CFastBTRemoveView)ON_COMMAND(ID_TEST1, OnTest1)ON_COMMAND(ID_WRITE, OnWrite)ON_COMMAND(ID_CORRECT, OnCorrect)ON_COMMAND(ID_ImageCorrect, OnImageCorrect)ON_COMMAND(ID_ALL, OnAll)ON_COMMAND(ID_FILE_OPEN, OnFileOpen)ON_WM_CREATE()ON_UPDATE_COMMAND_UI(ID_ALL, OnUpdateAll)/AFX_MSG_MA

22、P/ Standard printing commandsON_COMMAND(ID_FILE_PRINT, CView:OnFilePrint)ON_COMMAND(ID_FILE_PRINT_DIRECT, CView:OnFilePrint)ON_COMMAND(ID_FILE_PRINT_PREVIEW, CView:OnFilePrintPreview)END_MESSAGE_MAP()/ CFastBTRemoveView construction/destructionCFastBTRemoveView:CFastBTRemoveView()/ TODO: add constru

23、ction code hereCFastBTRemoveView:CFastBTRemoveView()BOOL CFastBTRemoveView:PreCreateWindow(CREATESTRUCT& cs)/ TODO: Modify the Window class or styles here by modifying/ the CREATESTRUCT csreturn CView:PreCreateWindow(cs);/ CFastBTRemoveView drawingvoid CFastBTRemoveView:OnDraw(CDC* pDC)CFastBTRemove

24、Doc* pDoc = GetDocument();ASSERT_VALID(pDoc);/ TODO: add draw code for native data herepDC-TextOut(0,Output_Y,Output_String);/ CFastBTRemoveView printingBOOL CFastBTRemoveView:OnPreparePrinting(CPrintInfo* pInfo)/ default preparationreturn DoPreparePrinting(pInfo);void CFastBTRemoveView:OnBeginPrint

25、ing(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)/ TODO: add extra initialization before printingvoid CFastBTRemoveView:OnEndPrinting(CDC* /*pDC*/, CPrintInfo* /*pInfo*/)/ TODO: add cleanup after printing/ CFastBTRemoveView diagnostics#ifdef _DEBUGvoid CFastBTRemoveView:AssertValid() constCView:AssertValid()

26、;void CFastBTRemoveView:Dump(CDumpContext& dc) constCView:Dump(dc);CFastBTRemoveDoc* CFastBTRemoveView:GetDocument() / non-debug version is inlineASSERT(m_pDocument-IsKindOf(RUNTIME_CLASS(CFastBTRemoveDoc);return (CFastBTRemoveDoc*)m_pDocument;#endif /_DEBUG/ CFastBTRemoveView message handlers/*数据读取

27、方式ret = theHDF.ReadQKMData(int band_index, int sx, int sy, int dx, int dy, int rgb)输入参数：band_index: 读取波段，对于QKM数据，只有BAND1Q,BAND2Qsx: 起始点左上角x坐标sy: 起始点左上角y坐标dx: 水平方向读取像元个数dy: 垂直方向读取像元个数rgb: 读入哪个通道，取值：1，2，3，分别代表R,G,B通道返回值：ret0: 数据读取成功-1: band_index参数错误-2: sx,sy,dx,dy参数为负数错误-3: rgb参数错误-4: sx(起始点x坐标)超过Np-

28、5: sy(起始点y坐标)超过Nl-7: SDreaddata出错R通道数据放在theHDF.pdataR (uint16*)G通道数据放在theHDF.pdataG (uint16*)B通道数据放在theHDF.pdataB (uint16*)取值范围为0-65535，其中0-32767是有效值*/void CFastBTRemoveView:OnTest1() / TODO: Add your command handler code hereCClientDC dc (this);dc.TextOut(0,0,Preparing for the data lodaing.);CTime

29、t; t = CTime:GetCurrentTime(); CString str = t.Format(%H:%M:%S); dc.TextOut(500,0,str);int ret;int retimage;ret = theHDF.ReadQKMData(LATQ, 0, 0, 1354, 4020, 0);/第一个0就是作为X坐标起点，第二个0不变。50是X范围，100是theHDF.img_height/4,最后的0保留。retimage = theHDF.ReadQKMData(BAND1Q, 0, 0, 5416, 16080, 1);void CFastBTRemoveVi

30、ew:OnWrite() / TODO: Add your command handler code hereCFile file;CString str;int i,j, k, index, index2, indeximage;/i代表每列循环，j代表每行循环,k代表第几个扫描行,index2代表内存临时数组记录顺序pLongLattemp = new float32 1354*4020;pImage = new uint16 5416*16080;CString teststr , testtemp ;CString *pOut;index2 = 0;indeximage = 0;for

31、(i=0; i4020; i+)for(j=0; j1354; j+)index = j + i * 1354; pLongLattempindex2+ = theHDF.pdataindex;/end of j/end of ifor (i=0;i16080;i+)for (j=0;j5416;j+)index=j+i*5416;pImageindeximage+ = theHDF.pdataRindex;/屏幕输出CClientDC dc (this);dc.TextOut(0,30,Data load completed.);CTime t; t = CTime:GetCurrentTi

32、me(); str = t.Format(%H:%M:%S); dc.TextOut(500,30,str); void CFastBTRemoveView:OnCorrect() / 去除蝴蝶结现象测试float32 old_data40;float32 new_data40;int num_old_data = 0; /0-3int num_new_data = 10; /0-9float32 interval; /new_data的间隔(0-3之间)float32 wrong_interval;float32 accu; /累积间隔float32 dist_to_prev, dist_t

33、o_next;int old_x;int i, j, k, m, n, index, a;int o, l; /用于输出检验变量int qver, qhor; /用于错误值的处理float32 temp;float32 tempup, tempdown;/用于错误值的处理CFile filecorr;CFile testfile;pLongLat_Num_Old_Data = new int 1354;for(j=10; j4020; j+=10)/垂直读取for(i=0; i=380 & i pLongLattemp index )temp = pLongLattempindex;for (

34、k=1; k=10; k+)if (pLongLattempindex - k * 1354 =0; l-)old_datal = pLongLattemp index - (l+k) * 1354;interval = (num_old_data-1)/(double)(num_new_data-1);new_data0 = old_data0; /新第0个=老第0个new_datanum_new_data-1 = old_datanum_old_data-1; /新最后个=老最后个old_x = 0;accu = 0;for(m=1; m=1)old_x+;accu = accu - 1;

35、 /总是0-1之间dist_to_prev = accu;dist_to_next = 1-accu;new_datam = old_dataold_x*dist_to_next + old_dataold_x+1*dist_to_prev;for(n=9; n=0; n-)pLongLattemp index-n*1354 = new_datan;break;/endfor i/endfor j/进行损坏数据的检验与纠正for (qver=0; qver4020; qver+)if (pLongLattempqver*1354-100)/开始对损坏数据进行纠正for (qhor=0;qhor

36、1354;qhor+)tempup=pLongLattempqver*1354+qhor-1354;tempdown=pLongLattempqver*1354+qhor+10*1354;wrong_interval=(tempdown-tempup)/11;for (l=0;l10;l+)pLongLattempqver*1354+qhor+l*1354=pLongLattempqver*1354+qhor+(l-1)*1354+wrong_interval; writef(E:纬度数据.hdf,纬度纠正数据,4020,1354,pLongLattemp);delete pLongLatte

37、mp;CClientDC dc (this);dc.TextOut(0,60,Latitude correction completed.);CTime t; t = CTime:GetCurrentTime(); CString str = t.Format(%H:%M:%S); dc.TextOut(500,60,str);dc.TextOut(0,90,Processing the image Bow-Tie correction, please WAIT.);t = CTime:GetCurrentTime(); str = t.Format(%H:%M:%S); dc.TextOut

38、(500,90,str);void CFastBTRemoveView:OnImageCorrect() / TODO: Add your command handler code hereuint16 old_data40;uint16 new_data40;int num_old_data = 0; /0-3int num_new_data = 40; /0-40float32 interval; /new_data的间隔(0-3之间)uint16 wrong_interval;float32 accu; /累积间隔float32 dist_to_prev, dist_to_next;in

39、t old_x;int i, j, k, m, n, index;int o, l,r; /用于输出检验变量int qver, qhor, num_wrong_data; /用于错误值的处理uint16 tempup, tempdown;/用于错误值的处理m_progBar.ShowWindow(TRUE);for(j=40; j4020*4; j+=40)/垂直读取m_progBar.StepIt();for(i=0; i=380*4 & i 0 )num_old_data = (pLongLat_Num_Old_Data(int)(i/4)*4;for(l=0; l=num_old_dat

40、a-1; l+)old_datal = pImage index - 40 * 1354 * 4 + l*1354*4; interval = (num_old_data-1)/(double)(num_new_data-1);new_data0 = old_data0; /新第0个=老第0个new_datanum_new_data-1 = old_datanum_old_data-1; /新最后个=老最后个 old_x = 0;accu = 0;for(m=1; m=1)old_x+;accu = accu - 1; /总是0-1之间dist_to_prev = accu;dist_to_n

41、ext = 1-accu;new_datam = old_dataold_x*dist_to_next + old_dataold_x+1*dist_to_prev;for(n=0; n=39; n+)pImage index - 40 * 1354 * 4 + n*1354*4 = new_datan;/end for i/end for j/ 进行损坏数据的检验与纠正for (qver=0; qver4020*4-40; qver+)m_progBar.StepIt();for (qhor=0;qhor=32767)/开始对损坏数据进行纠正for (l=0;l=32767)num_wron

42、g_data+;tempup=pImageqver*1354*4+ qhor - 1354*4;tempdown=pImageqver*1354*4+qhor+(num_wrong_data+1)*4*1354;wrong_interval=(tempdown-tempup)/(num_wrong_data+1);for (r=0;rnum_wrong_data+1;r+)pImageqver*1354*4+qhor+r*1354*4=pImageqver*4*1354+qhor+(r-1)*1354*4+wrong_interval;CClientDC dc (this);dc.TextOu

43、t(0,120,Writing HDF File.);CTime t; t = CTime:GetCurrentTime(); CString str = t.Format(%H:%M:%S); dc.TextOut(500,120,str);m_progBar.ShowWindow(FALSE); write(E:纠正图像.hdf,纠正图像,16080,5416,pImage);delete pImage;delete theHDF.pdata;delete theHDF.pdataR;dc.TextOut(0,150,Image data correction completed.);t

44、= CTime:GetCurrentTime(); str = t.Format(%H:%M:%S); dc.TextOut(500,150,str);void CFastBTRemoveView:write(char *filename, char *sdsname, int row, int col, void *pdataw)int32 sd_id, sds_id, sds_index;intn status;int32 start2, edges2;int32 dim_sizes2;/sds维度大小int i, j;/* end Variable declaration */ /* O

45、pen the file and initialize the SD interface.*/if(strcmp(sdsname,SDSred)=0)sd_id = SDstart (filename, DFACC_CREATE);else if(strcmp(sdsname,SDSgreen)=0|strcmp(sdsname,SDSblue)=0)sd_id = SDstart (filename, DFACC_WRITE);elsesd_id = SDstart (filename, DFACC_CREATE); CString s1;/s1.Format(sd_id is:%d,sd_

46、id);/MessageBox(s1);dim_sizes0=16080;dim_sizes1=5416;/* create a sds.*/ sds_id = SDcreate (sd_id, sdsname, DFNT_INT16, 2, dim_sizes);CString s;/s.Format(sds_id is:%d,sds_id);/MessageBox(s);/* * Define the location and size of the data to be written to the data set.*/start0 = 0;start1 = 0;edges0 = ro

47、w;edges1 = col;/* Write the stored data to the data set. The third argument is set to NULL* to specify contiguous data elements. The last argument must* be explicitly cast to a generic pointer since SDwritedata is designed* to write generic data. */if(status = SDwritedata (sds_id, start, NULL, edges

48、, pdataw)=-1) MessageBox(write data wrong!);exit(-1); /* Terminate access to the data set.*/status = SDendaccess (sds_id);/* Terminate access to the SD interface and close the file.*/status = SDend (sd_id);void CFastBTRemoveView:writef(char *filename, char *sdsname, int row, int col, float *pdataw)i

49、nt32 sd_id, sds_id, sds_index;intn status;int32 start2, edges2;int32 dim_sizes2;/sds维度大小int i, j;/* end Variable declaration */ /* Open the file and initialize the SD interface.*/if(strcmp(sdsname,SDSred)=0)sd_id = SDstart (filename, DFACC_CREATE);else if(strcmp(sdsname,SDSgreen)=0|strcmp(sdsname,SDSblue)=0)sd_id = SDstart (filename, DFACC_WRITE);elsesd_id = SDstart