畢業(yè)設計（論文）MODIS 1B影像的BowTie效應快速去除技術

1、modis 1b影像的bow-tie效應快速去除技術rapid removal of bow-tie effects from modis 1b data: an algorithmic development目錄摘要-3abstract-4一、 引言-5二、 算法設計-6（一） bow-tie現(xiàn)象成因分析-6（二） modis 1b的緯度科學數(shù)據(jù)集-7（三） bow-tie去除算法流程-8三、 算法實現(xiàn)-12（一） 緯度數(shù)據(jù)的蝴蝶結現(xiàn)象去除代碼示例-12（二） 圖像數(shù)據(jù)的蝴蝶結現(xiàn)象去除代碼示例-14四、 結果與分析-15五、 結論與討論-17參考文獻-19致謝-19附錄1程序運行界面-20附

2、錄2程序源代碼-21modis 1b影像的bow-tie效應快速去除技術摘 要：中分辨率成像光譜儀（modis）是一款廣泛應用于進行大、中尺度大氣、海洋和陸地綜合監(jiān)測的新一代被動式擺動掃描遙感傳感器，其55o的寬視場角和2330公里的幅寬使其影像邊緣區(qū)域存在明顯的重疊狀幾何畸變，稱為bow-tie效應或蝴蝶結效應。目前該類幾何畸變通常在進行投影轉換時同步去除，但投影轉換處理耗時較長，不利于影像的直接快速瀏覽。而現(xiàn)有bow-tie效應去除算法除了envi外多直接基于影像上重疊區(qū)的檢測，未充分利用modis自帶的經(jīng)緯度信息，算法效率不高。本文從modis影像的成像機制及其1b文件的存儲結構入手，根

3、據(jù)其自帶的星下點1公里間隔緯度數(shù)據(jù)計算掃描方向上各列的重疊程度，并以此為依據(jù)進行bow-tie效應去除。此外，鑒于部分緯度數(shù)據(jù)可能由于儀器等原因出現(xiàn)異常值，在本文的算法中也對這些異常進行了糾正。本算法在visual c+ 6.0環(huán)境下開發(fā)實現(xiàn)，數(shù)據(jù)處理速度與糾正后影像質量均達到了設計要求。關鍵詞：modis; bow-tie效應去除; 算法設計; 畸變?nèi)コ齬apid 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聯(lián)合美國和國際多個科研機構歷時8年、耗資13億美元研制的新一代星載遙感傳感器，其目標主要是建立全球大氣、海洋和陸地的動力學模型，并能在地球發(fā)生變化之前做出準確的預測。modis搭載在地球觀測系統(tǒng)(eos)的terra和aqua兩顆衛(wèi)星上 劉玉潔，楊忠東. modis遙感信息處理原理與算法m. 北京：科學出版社，2000.。eos計劃的目標是由一系列低軌道大型衛(wèi)星平臺裝載先進儀器，以獲得大量的先進衛(wèi)星遙感資料，進而通過反演獲得能反映地球這個復雜系統(tǒng)變異的多方面綜合信息。modis共有36個觀測通道，波

10、譜范圍從0.4m到14m，覆蓋范圍包括可見光、近紅外、短波紅外和遠紅外譜區(qū)。modis的星下點分辨率可分為250m、500m和1000m三種。目前modis共有44種標準產(chǎn)品，按研究對象可分為陸地、海洋和大氣三類；按接收和處理過程又可分為1a、1b、l2、l2g、l3和l4等多個級別 modis level 1a earth location: algorithm theoretical basis documente. 1997, 。其中經(jīng)輻射定標和地理坐標定位的modis 1b數(shù)據(jù)是本文的研究數(shù)據(jù)源。華東師范大學地理信息科學教育部重點實驗室

11、在自然地理館樓頂安裝了衛(wèi)星數(shù)字廣播接收系統(tǒng)(dvb-s)，可準實時接收國家衛(wèi)星氣象中心匯集處理后發(fā)布的modis 1b數(shù)據(jù)。本文所采用的示例數(shù)據(jù)為2006年10月9日國家衛(wèi)星氣象中心北京站接收的aqua/modis數(shù)據(jù)。目前在國內(nèi)外針對蝴蝶結現(xiàn)象去除的專門研究開展不多，絕大多數(shù)將蝴蝶結現(xiàn)象的去除與modis數(shù)據(jù)的幾何糾正同時進行。這樣雖然可直接利用現(xiàn)有幾何糾正算法，如多項式糾正模型等，但是其執(zhí)行效率卻非常低。國外針對modis 1b數(shù)據(jù)的處理軟件mrt(free modis reprojection tool)以及通用遙感影像處理軟件envi均可以進行這樣的處理（注：envi 4.2版本后中的

12、相關算法已開始基于地理坐標進行處理，但仍未將蝴蝶結效應去除與圖像及和糾正獨立出來），但根據(jù)處理的數(shù)據(jù)分辨率不同，其耗費的時間往往需要數(shù)分鐘甚至是幾十分鐘。對于單獨的蝴蝶結現(xiàn)象去除的研究相對較少。郭廣猛(2003)提出了利用idl語言開辟一個33的滑動窗口，來對整幅影像進行掃描判斷改矩陣的相關系數(shù)，并以此對重疊的圖像進行糾正；程亮等(2005)利用modis 1b數(shù)據(jù)所對應的mod03數(shù)據(jù)即地理修訂文件來對圖像進行掃描，并利用分形插值方法進行蝴蝶結現(xiàn)象的去除。其他的相關研究基本也是圍繞著兩種方法展開。這兩種方法或者僅從圖像學的角度進行蝴蝶結現(xiàn)象的去除或者利用modis的星歷文件進行去除，算法設計

13、較為復雜，執(zhí)行效率也相對較低。因此在modis數(shù)據(jù)接收終端等應用環(huán)境下，如何快速地對蝴蝶結現(xiàn)象進行去除，增強影像的可讀性便成為本文的主要研究目的。二、 算法設計（一） bow-tie現(xiàn)象成因分析modis探測器是一種被動式擺動掃描探測器，其橫向掃描角為55。每條掃描線由1 354個1 km的像素或2 708個500 m像素或5 416個250 m像素組成。由于地球曲率的影響,掃描線的實際跨度大約為2 330km。每完成一次掃描,modis探測器沿軌道前進10 km,這10 km的航向區(qū)域就構成了一個掃描條帶。每個掃描條帶的寬度分別為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數(shù)據(jù)存在多種幾何畸變, 其中modis 1b數(shù)據(jù)相鄰掃描帶之間存在的錯位和重疊現(xiàn)象（即蝴蝶結效應）是其重要的系統(tǒng)畸變之一。“蝴蝶結”

15、現(xiàn)象隨著觀測角度的增大而趨嚴重,計算表明,在視場角為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圖像重疊現(xiàn)象的方法研究j. 山西氣象, 2005, (3), p1-3。圖1modis 1b 1km數(shù)據(jù)左半部分蝴蝶結效應示意圖figure 1. the bow-tie effect on left part of the modis 1b 1km data通過上

17、述分析可見，正是由于modis的這種成像機制，決定了其影像不可避免地帶有蝴蝶結效應，且該效應越往兩側越嚴重，這種系統(tǒng)性幾何畸變大大影響了影像的可讀性，尤其是影像邊緣區(qū)的重疊現(xiàn)象幾乎完全掩蓋了地表的真實情況，對modis 1b數(shù)據(jù)蝴蝶結效應的去除是對modis數(shù)據(jù)進行預處理的必備操作。（二） modis 1b的科學數(shù)據(jù)集modis 1b數(shù)據(jù)主要以hdf文件的格式進行存儲和分發(fā)。hdf（hierarchical data format）格式是由美國伊利諾大學的ncsa（national center for supercomputing applications）組織開發(fā)的一種超文本文件格式，是一

18、種用于存儲和分發(fā)科學數(shù)據(jù)的自描述性、多對象的層次數(shù)據(jù)格式 景毅剛 等，hdf-eos數(shù)據(jù)格式及其數(shù)據(jù)導入j, 氣象科技，2005-6 333，p278-282。通常其記錄的對象包含多維排列,用于存儲科學數(shù)據(jù)，即稱為科學數(shù)據(jù)集（sds）。1b級使用sds存儲校準過的科學數(shù)據(jù),它們的不確定性,部分質量保證數(shù)據(jù)。通過屬性數(shù)據(jù)使sdss能夠自我描述。利用hdf所提供的相關庫函數(shù)可以很方便將文件中存儲的原始數(shù)據(jù)讀入內(nèi)存。由于modis 1b數(shù)據(jù)自帶有經(jīng)緯度信息，并以1km分辨率存于1km、500m、250m各個分辨率的數(shù)據(jù)文件中 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/*數(shù)據(jù)讀取

27、方式ret = thehdf.readqkmdata(int band_index, int sx, int sy, int dx, int dy, int rgb)輸入?yún)?shù)：band_index: 讀取波段，對于qkm數(shù)據(jù)，只有band1q,band2qsx: 起始點左上角x坐標sy: 起始點左上角y坐標dx: 水平方向讀取像元個數(shù)dy: 垂直方向讀取像元個數(shù)rgb: 讀入哪個通道，取值：1，2，3，分別代表r,g,b通道返回值：ret0: 數(shù)據(jù)讀取成功-1: band_index參數(shù)錯誤-2: sx,sy,dx,dy參數(shù)為負數(shù)錯誤-3: rgb參數(shù)錯誤-4: sx(起始點x坐標)超過np-

28、5: sy(起始點y坐標)超過nl-7: sdreaddata出錯r通道數(shù)據(jù)放在thehdf.pdatar (uint16*)g通道數(shù)據(jù)放在thehdf.pdatag (uint16*)b通道數(shù)據(jù)放在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代表每列循環(huán)，j代表每行循環(huán),k代表第幾個掃描行,index2代表內(nèi)存臨時數(shù)組記錄順序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() / 去除蝴蝶結現(xiàn)象測試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/進行損壞數(shù)據(jù)的檢驗與糾正for (qver=0; qver4020; qver+)if (plonglattempqver*1354-100)/開始對損壞數(shù)據(jù)進行糾正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:緯度數(shù)據(jù).hdf,緯度糾正數(shù)據(jù),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/ 進行損壞數(shù)據(jù)的檢驗與糾正for (qver=0; qver4020*4-40; qver+)m_progbar.stepit();for (qhor=0;qhor=32767)/開始對損壞數(shù)據(jù)進行糾正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 */ /* open the file and initialize the sd interface.*/if(strcmp(sdsname,sdsred)=0)sd_id = sdstart (