字全息離軸干涉:從原理到代碼實現(xiàn))
簡介本資源是一套面向光學(xué)工程、數(shù)字全息與計算成像方向的本科高年級或研究生實踐教學(xué)材料聚焦像面數(shù)字全息中離軸干涉圖樣的仿真建模與相位恢復(fù)全流程。資源完整實現(xiàn)從物光波前生成、參考光引入、干涉圖模擬、頻譜濾波、傅里葉逆變換到相位解包裹與強度再現(xiàn)的全部核心算法適用于課程設(shè)計、實驗復(fù)現(xiàn)及算法原理驗證。壓縮包共24個文件含6個Matlab源碼如fresnel.m、unwrap_LS.m、New_FFT2.m等關(guān)鍵函數(shù)、13張運行結(jié)果圖jpg格式涵蓋原始干涉圖、頻譜圖、相位像、強度再現(xiàn)像等中間與最終可視化結(jié)果以及5幅bmp格式的原始/處理圖像如usaf.bmp、0相位像.bmp整體僅869KB輕量易用。已有476人學(xué)習(xí)下載配套代碼結(jié)構(gòu)清晰、變量命名規(guī)范、關(guān)鍵步驟附注釋可直接運行并支持參數(shù)調(diào)整特別適合作為數(shù)字全息入門算法理解與Matlab光學(xué)仿真能力提升的實操范例。1. 項目概述從一束光到一幅圖如果你對光學(xué)、圖像處理或者計算成像感興趣手里又恰好有一份名為“模擬像面數(shù)字全息離軸干涉含Matlab源碼.zip”的文件那你算是來對地方了。這不僅僅是一個壓縮包它更像是一把鑰匙能幫你打開一扇通往“計算光學(xué)”這個交叉領(lǐng)域的大門。簡單來說這個項目就是用計算機軟件Matlab來模擬和重現(xiàn)一個經(jīng)典的光學(xué)實驗過程——離軸數(shù)字全息干涉并且是在“像面”這個特殊的位置上進(jìn)行。聽起來有點繞別急我干了十幾年光學(xué)仿真和圖像處理今天就用最“人話”的方式帶你把這個項目從里到外、從原理到代碼徹底盤明白。數(shù)字全息你可以把它理解成給光波“拍一張包含所有信息的立體照片”。傳統(tǒng)相機只能記錄光的強度也就是明暗但全息技術(shù)能同時記錄下光的強度和相位你可以理解為光波的“形狀”或“起伏”。離軸干涉則是實現(xiàn)這種記錄的一種巧妙方法它通過引入一個傾斜的參考光把物體的相位信息“編碼”到一幅我們能直接拍攝到的干涉條紋圖中。而“模擬像面”指的是我們模擬的這個全息圖不是在物體所在的平面記錄的而是在經(jīng)過一個透鏡成像后的像平面上記錄的。這么做的好處是物體的像和全息圖在同一平面上重建過程更直觀也更貼近很多實際顯微成像的應(yīng)用場景。這個Matlab源碼包的價值在于它把一套復(fù)雜的光學(xué)理論和數(shù)學(xué)運算封裝成了你可以一行行運行、一步步修改的代碼。你不需要昂貴的激光器、精密的光學(xué)平臺和防震臺只需要一臺電腦就能親手“搭建”這個光學(xué)系統(tǒng)觀察干涉圖樣的生成并最終通過計算“再現(xiàn)”出物體的三維光場。無論是學(xué)生做課程設(shè)計、研究者驗證算法還是工程師進(jìn)行系統(tǒng)前期仿真這都是一個極佳的起點。接下來我會帶你拆解它的核心思路、手把手過一遍代碼、分享我踩過的坑和調(diào)試技巧讓你不僅能運行出結(jié)果更能理解每一步背后的“所以然”。2. 核心原理拆解光是如何被“計算”的要玩轉(zhuǎn)這個源碼我們不能只當(dāng)個“調(diào)包俠”或者“跑碼員”必須得弄清楚它到底在算些什么。這里面的核心物理是“干涉”核心數(shù)學(xué)是“傅里葉變換”和“角譜傳播”核心操作是“模擬”和“重建”。我們一層層來看。2.1 離軸全息干涉的編碼原理為什么離軸干涉這么重要在早期的同軸全息參考光和物光方向一致中重建的物體像會和共軛像、零級衍射光重疊在一起畫面非?;靵y。離軸技術(shù)就像在合唱中讓領(lǐng)唱站到一邊他的聲音物光信息就能從和聲零級光和其他噪聲中被清晰地分辨出來。在模擬中我們通常用兩個復(fù)數(shù)矩陣來表示物光波 ( O(x, y) ) 和參考光波 ( R(x, y) )。物光波攜帶著物體信息它的振幅和相位由物體決定。參考光波通常被模擬成一個理想的平面波其數(shù)學(xué)表達(dá)式為 ( R(x, y) A_r \exp(i 2\pi (f_x x f_y y)) )。這里的 ( A_r ) 是常數(shù)振幅關(guān)鍵就在于 ( (f_x, f_y) )它們是參考光的空間頻率決定了光波的傾斜角度。正是這個傾斜使得物光和參考光干涉后物體的相位信息被轉(zhuǎn)換成了干涉條紋的疏密和走向。兩者在像面上疊加其光強分布也就是我們模擬能“拍”到的全息圖 ( I_H(x, y) )為 [ I_H |O R|^2 |O|^2 |R|^2 O^R OR^] 這個公式是全息的靈魂。( |O|^2 ) 是物光自身強度( |R|^2 ) 是參考光強度它們合起來構(gòu)成緩慢變化的背景零級項。而 ( O^R ) 和 ( OR^) 是交叉項其中就編碼了我們需要的物體相位信息。在離軸設(shè)置下這兩個交叉項在頻域傅里葉變換后的空間里會和零級項分離開這是后續(xù)我們能干凈利落地提取出物體信息的前提。2.2 像面全息與角譜傳播為什么強調(diào)“像面”在菲涅耳全息或夫瑯禾費全息中全息圖記錄面距離物體較遠(yuǎn)。而在像面全息中我們假設(shè)物體已經(jīng)通過一個理想透鏡完美成像在了記錄平面上。這意味著記錄平面上的物光復(fù)振幅 ( O(x, y) )直接就是物體透過率函數(shù)經(jīng)過透鏡成像后的結(jié)果無需考慮復(fù)雜的衍射傳播過程。這極大地簡化了模擬的物理模型。但是當(dāng)我們從全息圖重建物體光場時我們通常需要“回溯”光場到物體平面或者至少到某個觀察平面。這時就需要用到“角譜傳播”理論。角譜傳播法是計算光波在自由空間中傳播的精確方法之一在標(biāo)量衍射理論范圍內(nèi)。其核心思想是任何一個光波場都可以分解成無數(shù)個不同方向傳播的平面波角譜。傳播一段距離 ( z ) 后每個平面波成分只會產(chǎn)生一個相位延遲。在頻域里這個操作變得異常簡單。設(shè)輸入平面光場為 ( U(x, y; 0) )其二維傅里葉變換為 ( \mathcal{F}{U(x, y; 0)} G(f_x, f_y; 0) )其中 ( f_x, f_y ) 是空間頻率。那么傳播距離 ( z ) 后的光場頻譜為 [ G(f_x, f_y; z) G(f_x, f_y; 0) \cdot H(f_x, f_y) ] 其中傳遞函數(shù) ( H(f_x, f_y) \exp(i k z \sqrt{1 - (\lambda f_x)^2 - (\lambda f_y)^2}) )這里 ( k2\pi/\lambda ) 是波數(shù)( \lambda ) 是光波長。這個公式是源碼中重建部分的核心。通過一次傅里葉變換、乘以一個相位因子、再做一次逆傅里葉變換我們就完成了光場的傳播模擬。這種方法特別適合在計算機上實現(xiàn)也是數(shù)字全息重建的主流算法。2.3 數(shù)字重建的四步流程源碼的算法骨架無論是模擬還是重建都遵循一個清晰的四步流程理解了它看代碼就事半功倍模擬生成全息圖構(gòu)建或加載物體圖像 - 賦予其相位分布模擬物體形貌- 將其作為物光波 - 模擬透鏡成像效應(yīng)可能包含衍射- 與傾斜的參考光干涉 - 計算光強得到模擬全息圖。這一步的輸出是一幅灰度圖像。頻域濾波與信息提取對模擬全息圖進(jìn)行二維傅里葉變換得到其頻譜。在頻譜圖中你會看到中心亮團(tuán)零級項和兩側(cè)的“邊帶”交叉項。通過設(shè)計一個合適的濾波器通常是圓形或矩形窗函數(shù)將其中一個邊帶比如 ( O^*R ) 項完整地截取出來并移動到頻譜中心。這一步的目的是分離出攜帶物體信息的單項。角譜傳播重建對濾波并移中后的頻譜進(jìn)行逆傅里葉變換得到復(fù)數(shù)形式的像面光場分布。如果需要觀察不同深度的重建像就利用上面提到的角譜傳播公式將這個光場向不同距離 ( z ) 進(jìn)行“回溯”或“聚焦”計算。結(jié)果顯示與解包裹計算重建光場的振幅和相位。振幅圖就是物體的強度像比較直觀。相位圖則包含了物體的表面高度或折射率變化信息。但直接計算出的相位值被包裹在 ( [-\pi, \pi] ) 的主值區(qū)間內(nèi)呈鋸齒狀需要通過“相位解包裹”算法將其恢復(fù)成連續(xù)的相位分布這才是最終可用的三維形貌數(shù)據(jù)。注意在像面全息模擬中第一步的“透鏡成像效應(yīng)”有時會被簡化甚至直接認(rèn)為物光波就是物體本身如果物體是二維圖像。而在重建時如果只關(guān)心像面本身傳播距離 ( z ) 可以設(shè)為0。源碼的靈活性往往就體現(xiàn)在這些地方。3. 源碼結(jié)構(gòu)深度解析與關(guān)鍵函數(shù)剖析拿到一個“含Matlab源碼”的壓縮包最忌諱的就是直接運行main.m然后看個結(jié)果圖了事。我們必須像解剖一樣看清它的每一塊肌肉和骨骼。下面我以一個典型的、結(jié)構(gòu)良好的模擬像面離軸全息Matlab項目為例拆解其核心文件和作用。你的源碼包可能文件名不同但功能模塊大抵如此。3.1 主腳本 (main.m或demo_simulate_reconstruct.m)這是程序的入口通常負(fù)責(zé)設(shè)置全局參數(shù)、調(diào)用各個功能函數(shù)、控制流程和顯示結(jié)果。一個好的主腳本應(yīng)該像一份實驗報告提綱。%% 1. 清空與初始化 clear; close all; clc; addpath(genpath(./utils)); % 添加工具函數(shù)路徑 %% 2. 系統(tǒng)參數(shù)設(shè)置 lambda 632.8e-9; % 波長單位米 (He-Ne激光) k 2*pi/lambda; % 波數(shù) pixel_size 6.5e-6; % CCD像元尺寸單位米 Nx 1024; Ny 1024; % 全息圖尺寸像素 Lx Nx * pixel_size; % 全息圖物理寬度 Ly Ny * pixel_size; % 全息圖物理高度 % 參考光傾斜角度通常用離軸頻率表示 theta_x 1.0 * pi/180; % x方向傾斜角弧度 theta_y 0.5 * pi/180; % y方向傾斜角弧度 fx_r sin(theta_x) / lambda; % 參考光空間頻率x分量 fy_r sin(theta_y) / lambda; % 參考光空間頻率y分量 %% 3. 生成或加載物體 object_amp im2double(imread(test_object.png)); % 物體振幅 object_phase generate_phase_map(Nx, Ny); % 物體相位模擬形貌 object_field object_amp .* exp(1i * object_phase); % 物光復(fù)振幅 %% 4. 模擬像面全息圖記錄 [hologram, ref_field] simulate_imgholo(object_field, pixel_size, lambda, fx_r, fy_r); %% 5. 全息圖重建 recon_distance 0; % 像面重建距離為0 [recon_amp, recon_phase] reconstruct_imgholo(hologram, pixel_size, lambda, fx_r, fy_r, recon_distance); %% 6. 結(jié)果顯示 figure; subplot(2,3,1); imshow(object_amp, []); title(原始物體振幅); subplot(2,3,2); imshow(object_phase, []); title(原始物體相位); subplot(2,3,3); imshow(hologram, []); title(模擬像面全息圖); subplot(2,3,4); imshow(abs(fftshift(fft2(hologram))), []); title(全息圖頻譜對數(shù)); subplot(2,3,5); imshow(recon_amp, []); title(重建振幅); subplot(2,3,6); imshow(recon_phase, []); title(重建包裹相位);在主腳本中參數(shù)設(shè)置部分是重中之重。pixel_size和lambda的比例決定了系統(tǒng)的最大記錄角度。fx_r,fy_r的設(shè)置必須足夠大以確保頻譜中的各項充分分離但又不能太大否則干涉條紋過密超過CCD的奈奎斯特采樣頻率會導(dǎo)致混疊。一個經(jīng)驗法則是參考光引入的載頻應(yīng)使邊帶中心與頻譜中心的距離大于物體帶寬的3倍左右。3.2 全息圖模擬函數(shù) (simulate_imgholo.m)這個函數(shù)是物理模型的核心。它精確地模擬了物光傳播到像面、與參考光干涉的過程。function [hologram, ref_field] simulate_imgholo(object_field, pixel_size, lambda, fx_r, fy_r) % 輸入物光場像素尺寸波長參考光空間頻率 % 輸出全息圖強度參考光場可選 [Ny, Nx] size(object_field); Lx Nx * pixel_size; Ly Ny * pixel_size; % 生成坐標(biāo)網(wǎng)格 x linspace(-Lx/2, Lx/2, Nx); y linspace(-Ly/2, Ly/2, Ny); [X, Y] meshgrid(x, y); % 1. 模擬像面物光場這里簡化處理認(rèn)為物體已完美成像在像面。 % 更復(fù)雜的模擬可以在此處加入透鏡傳遞函數(shù)或衍射傳播。 img_field object_field; % 2. 生成傾斜平面波作為參考光 % 關(guān)鍵參考光相位是線性相移對應(yīng)空間頻率 ref_field exp(1i * 2*pi * (fx_r * X fy_r * Y)); % 3. 干涉并記錄強度 interference_field img_field ref_field; hologram abs(interference_field).^2; % 4. (可選) 添加噪聲更接近真實實驗 % hologram imnoise(hologram, gaussian, 0, 0.01); end關(guān)鍵點剖析meshgrid生成的X,Y是物理坐標(biāo)網(wǎng)格單位是米。這是所有物理計算的基礎(chǔ)千萬不能和像素索引i, j混淆。ref_field的生成是離軸干涉的靈魂代碼。fx_r * X fy_r * Y計算了每個坐標(biāo)點相對于原點的光程差乘以2π后就是相位延遲。這個線性相位項在傅里葉變換后會表現(xiàn)為頻譜的平移。干涉強度hologram是一個實數(shù)矩陣取值范圍從0到最大值。在顯示時通常需要歸一化或調(diào)整對比度。3.3 全息圖重建函數(shù) (reconstruct_imgholo.m)這是數(shù)字全息的魔法所在將一幅干涉條紋圖還原為三維光場。function [recon_amp, recon_phase] reconstruct_imgholo(hologram, pixel_size, lambda, fx_r, fy_r, z) % 輸入全息圖像素尺寸波長參考光頻重建距離z % 輸出重建平面的振幅和相位 [Ny, Nx] size(hologram); Lx Nx * pixel_size; Ly Ny * pixel_size; % 1. 全息圖頻譜分析 H fft2(hologram); H_shift fftshift(H); % 將零頻移到中心方便觀察 % 2. 設(shè)計頻域濾波器提取1級邊帶 % 生成頻域坐標(biāo) fx linspace(-1/(2*pixel_size), 1/(2*pixel_size), Nx); fy linspace(-1/(2*pixel_size), 1/(2*pixel_size), Ny); [Fx, Fy] meshgrid(fx, fy); % 預(yù)測1級邊帶中心位置對應(yīng)O*R項 fx_center -fx_r; % 注意符號由干涉公式推導(dǎo)決定 fy_center -fy_r; % 創(chuàng)建圓形濾波器也可用矩形窗 radius min(1/(4*pixel_size), sqrt(fx_center^2 fy_center^2)/3); % 濾波器半徑經(jīng)驗值 filter_mask ((Fx - fx_center).^2 (Fy - fy_center).^2) radius^2; % 3. 濾波并移中 H_filtered H_shift .* filter_mask; % 將濾波后的邊帶移回頻譜原點以消除參考光傾斜引入的線性相位 H_centered ifftshift(H_filtered); % 4. 逆傅里葉變換得到復(fù)振幅 recon_field_complex ifft2(H_centered); % 5. 角譜傳播如果z不為0 if z ~ 0 % 生成角譜傳遞函數(shù)的頻率坐標(biāo)注意與濾波時頻率坐標(biāo)的尺度可能不同 dfx 1 / Lx; dfy 1 / Ly; u ((0:Nx-1) - floor(Nx/2)) * dfx; % 注意此處的頻率坐標(biāo)排列需與fft2匹配 v ((0:Ny-1) - floor(Ny/2)) * dfy; [U, V] meshgrid(u, v); % 角譜傳遞函數(shù) H_as exp(1i * 2*pi/lambda * z .* sqrt(1 - (lambda*U).^2 - (lambda*V).^2)); % 傳播操作 recon_field_complex ifft2( fft2(recon_field_complex) .* ifftshift(H_as) ); end % 6. 提取振幅和相位 recon_amp abs(recon_field_complex); recon_phase angle(recon_field_complex); % 這是包裹相位 end這是最容易出錯的模塊需要逐行理解頻譜坐標(biāo)生成linspace(-1/(2*pixel_size), 1/(2*pixel_size), Nx)生成的是物理頻率坐標(biāo)單位是m^{-1}。其范圍由采樣定理決定最大頻率奈奎斯特頻率是1/(2*pixel_size)。這是設(shè)計濾波器的基準(zhǔn)。邊帶中心預(yù)測fx_center -fx_r這個符號關(guān)系需要根據(jù)具體干涉公式推導(dǎo)來驗證。一個快速驗證方法是運行模擬后查看全息圖的頻譜圖 (imshow(log(abs(H_shift)1), []))觀察亮斑位置然后調(diào)整符號使濾波器能對準(zhǔn)它。濾波器設(shè)計圓形濾波器的半徑radius是關(guān)鍵參數(shù)。太小會丟失物體高頻信息重建像模糊太大會包含零級或其他邊帶成分引入噪聲。通常需要根據(jù)頻譜圖手動調(diào)整。高級方法可以用自適應(yīng)濾波。角譜傳播的頻率坐標(biāo)注意第5步中生成U, V的方式與第2步生成Fx, Fy的方式可能不同。這是因為fft2輸出的頻率順序是[0, 1, ..., N/2, -N/21, ..., -1]。使用ifftshift(H_as)是為了將傳遞函數(shù)的零頻調(diào)整到與fft2輸出匹配的位置。這一步是角譜傳播正確實現(xiàn)的核心細(xì)節(jié)很多開源代碼在這里都有細(xì)微差別容易導(dǎo)致重建像散焦或出現(xiàn)偽影。3.4 工具函數(shù)集 (utils/文件夾)一個完整的項目通常還包含一些工具函數(shù)例如generate_phase_map.m: 生成模擬的物體相位分布如球形相位模擬微透鏡、斜坡相位模擬傾斜平面或隨機相位板。phase_unwrapping.m: 實現(xiàn)相位解包裹算法如 Goldstein枝切法、最小二乘法或基于質(zhì)量圖導(dǎo)向的路徑積分法。remove_tilt.m: 從重建相位中移除由于系統(tǒng)未對準(zhǔn)可能殘留的線性相位傾斜。calculate_metrics.m: 計算重建圖像與原始圖像的誤差指標(biāo)如均方誤差MSE、結(jié)構(gòu)相似性SSIM等。4. 從零開始手把手實現(xiàn)與參數(shù)調(diào)試實戰(zhàn)理解了原理和代碼結(jié)構(gòu)我們現(xiàn)在可以自己動手從一個簡單的例子開始逐步增加復(fù)雜度并學(xué)會如何調(diào)試。4.1 第一步創(chuàng)建一個簡單的二維振幅物體我們從一個最簡單的“雙縫”或“USAF分辨率板”圖像開始。這能讓我們清晰地看到干涉和重建效果。% 生成一個簡單的振幅物體 N 512; object_amp zeros(N, N); object_amp(200:220, 100:400) 1; % 一個水平亮條 object_amp(300:320, 100:400) 1; % 另一個水平亮條 % 或者加載一個標(biāo)準(zhǔn)分辨率板圖像 % object_amp im2double(imread(usaf1951.png)); % object_amp imresize(object_amp, [N, N]); % 初始相位設(shè)為0純振幅物體 object_phase zeros(N, N); object_field object_amp .* exp(1i * object_phase);運行主腳本你應(yīng)該能看到清晰的干涉條紋全息圖其頻譜呈現(xiàn)對稱的三個亮斑中心零級和兩側(cè)邊帶。重建的振幅圖應(yīng)該能大致恢復(fù)出兩條亮帶。常見問題1全息圖沒有條紋檢查參考光空間頻率fx_r,fy_r是否設(shè)置得太小嘗試增大角度例如從0.1度增加到1度。檢查物光場object_field是否為復(fù)數(shù)純實數(shù)只有振幅的物光與參考光干涉條紋對比度可能很低。確保有相位項exp(1i*phase)。4.2 第二步引入相位物體并觀察頻譜分離現(xiàn)在我們給物體加上相位信息模擬一個更真實的情況比如一個凸起的圓形相位物體。% 生成一個圓形相位物體 [X, Y] meshgrid(1:N, 1:N); center N/2; radius 50; phase_object zeros(N, N); ind sqrt((X-center).^2 (Y-center).^2) radius; phase_object(ind) 2 * pi * 2; % 給圓形區(qū)域一個2*2π的相位延遲一個“相位丘” object_field object_amp .* exp(1i * phase_object); % 結(jié)合振幅和相位再次運行。觀察全息圖的頻譜 (fftshift(fft2(hologram)))。你會發(fā)現(xiàn)兩側(cè)的邊帶亮斑不再是干凈的一個點而是擴(kuò)散成了一個“云團(tuán)”。這個云團(tuán)的大小就代表了物體信息的帶寬。確保你的濾波器半徑大于這個云團(tuán)的半徑否則會丟失高頻細(xì)節(jié)重建的相位圖邊緣會變得平滑。4.3 第三步精細(xì)調(diào)整濾波器與重建距離這是獲得高質(zhì)量重建結(jié)果的關(guān)鍵步驟也是最需要耐心和經(jīng)驗的地方??梢暬l譜并手動選取濾波器H_spectrum log(abs(fftshift(fft2(hologram))) 1); figure; imshow(H_spectrum, []); title(全息圖對數(shù)頻譜); % 交互式選取濾波器中心和大小的代碼簡化示意 % [fx_center, fy_center] ginput(1); % 用鼠標(biāo)點擊邊帶中心 % radius input(請輸入濾波器半徑像素: );在實際調(diào)試中我經(jīng)常先用ginput函數(shù)手動在頻譜圖上點選邊帶中心獲取其像素坐標(biāo)再換算成物理頻率坐標(biāo)。這樣可以精確校準(zhǔn)fx_center,fy_center。優(yōu)化濾波器形狀圓形濾波器簡單但可能會切掉方形成分。可以嘗試高斯濾波器或巴特沃斯濾波器讓濾波邊緣更平滑減少吉布斯振鈴現(xiàn)象。% 高斯濾波器示例 sigma radius / 2; % 高斯標(biāo)準(zhǔn)差 gaussian_mask exp(-((Fx - fx_center).^2 (Fy - fy_center).^2) / (2*sigma^2)); filter_mask gaussian_mask;調(diào)整重建距離z對于像面全息z0通常就是最清晰的平面。但你可以嘗試微調(diào)z比如正負(fù)幾毫米觀察重建振幅圖像的清晰度變化。這可以幫你理解“數(shù)字聚焦”的過程。編寫一個循環(huán)遍歷一系列z值計算每個z下重建圖像的銳度例如梯度平方和自動找到最佳聚焦位置。4.4 第四步相位解包裹與結(jié)果分析直接angle()函數(shù)得到的相位是包裹的對于我們的圓形相位丘你會看到一個從-π到π反復(fù)跳變的同心圓環(huán)。% 使用簡單的二維相位解包裹算法例如基于質(zhì)量圖導(dǎo)向的路徑積分法 % 假設(shè)有 phase_unwrapping 函數(shù) unwrapped_phase phase_unwrapping(recon_phase); % 可視化對比 figure; subplot(1,2,1); imshow(recon_phase, [-pi pi]); title(包裹相位); colorbar; subplot(1,2,2); imshow(unwrapped_phase, []); title(解包裹相位); colorbar; % 繪制中心線剖面 center_line unwrapped_phase(N/2, :); figure; plot(center_line); title(解包裹相位中心線剖面); xlabel(像素); ylabel(相位弧度);解包裹后的相位圖應(yīng)該呈現(xiàn)一個平滑的“丘”狀。你可以測量丘的高度相位值與我們最初設(shè)定的2 * 2π進(jìn)行對比驗證重建的保真度。如果存在整體傾斜或彎曲可能是系統(tǒng)像差或參考光不完全準(zhǔn)直導(dǎo)致的可以使用remove_tilt函數(shù)進(jìn)行校正。5. 避坑指南與性能優(yōu)化來自一線的經(jīng)驗在實際操作和代碼優(yōu)化中我積累了一些教科書和論文里很少提及但至關(guān)重要的經(jīng)驗。5.1 頻譜混疊與采樣定理這是數(shù)字全息模擬中最常見的“坑”。現(xiàn)象是重建圖像出現(xiàn)奇怪的波紋、扭曲或者高頻部分完全丟失。根源參考光傾斜角太大導(dǎo)致干涉條紋的空間頻率超過了奈奎斯特極限即每像素少于2個采樣點。診斷檢查全息圖本身。用plot(hologram(512, :))畫一行剖面線數(shù)一數(shù)條紋的周期。最小周期兩個亮紋或暗紋之間的距離對應(yīng)的像素數(shù)應(yīng)大于2。解決降低參考光角度減小fx_r,fy_r。增大像素尺寸在模擬中這意味著你的CCD像元更大雖然不現(xiàn)實但可以作為驗證。增加插值在模擬生成物光場和參考光場時使用更高的采樣率比如2倍生成全息圖后再降采樣回來。這相當(dāng)于過采樣。5.2 能量均衡與動態(tài)范圍問題全息圖頻譜中零級項中心亮斑的能量遠(yuǎn)高于邊帶導(dǎo)致濾波時邊帶信號微弱重建圖像信噪比低。解決調(diào)整物參光強比在模擬中可以給物光場乘以一個系數(shù)例如object_field 0.3 * object_field使物光與參考光強度匹配。經(jīng)驗上物參光強比在1:1到1:10之間效果較好。頻譜濾波前進(jìn)行背景抑制從全息圖中減去其平均強度或進(jìn)行高通濾波可以一定程度上壓制零級。5.3 計算效率與內(nèi)存優(yōu)化當(dāng)圖像尺寸很大如2048x2048或需要批量處理時計算速度和內(nèi)存成為瓶頸。使用單精度Matlab默認(rèn)雙精度。對于大部分全息計算單精度浮點數(shù) (single) 足夠且能節(jié)省一半內(nèi)存提升計算速度。在數(shù)據(jù)初始化時使用single()。object_field single(object_amp) .* exp(1i * single(object_phase));使用GPU加速如果你的Matlab支持Parallel Computing Toolbox且有好顯卡將大數(shù)據(jù)轉(zhuǎn)移到GPU上計算是質(zhì)的飛躍。object_field_gpu gpuArray(object_field); hologram_gpu abs(object_field_gpu ref_field_gpu).^2; hologram gather(hologram_gpu);注意頻繁在CPU和GPU之間傳輸數(shù)據(jù)會成為新的瓶頸盡量讓整個計算流程都在GPU上完成。優(yōu)化循環(huán)與向量化盡量避免在大的for循環(huán)中進(jìn)行逐像素操作。Matlab的矩陣運算向量化要快得多。我們之前的代碼已經(jīng)完全是向量化的。5.4 與實驗數(shù)據(jù)的對接這個模擬項目的終極價值之一是指導(dǎo)真實實驗。當(dāng)你用這套代碼處理實驗采集的全息圖時要注意像素尺寸校準(zhǔn)實驗中的pixel_size必須精確測量相機像元尺寸除以物鏡放大倍數(shù)。一個像素的誤差會導(dǎo)致重建尺度錯誤。參考光頻率標(biāo)定實驗中的fx_r,fy_r很難直接測量。最實用的方法是從實驗全息圖的頻譜中直接讀取。用鼠標(biāo)點選頻譜圖中邊帶亮斑的中心像素坐標(biāo)(cx, cy)然后通過公式換算fx_r_exp -(cx - floor(Nx/2) - 1) / (Nx * pixel_size); % 注意坐標(biāo)偏移和符號 fy_r_exp -(cy - floor(Ny/2) - 1) / (Ny * pixel_size);消除零級和共軛像實驗全息圖的頻譜分離可能不如模擬的完美。除了濾波還可以采用相位偏移法Phase Shifting等記錄多幅全息圖來徹底消除零級和共軛像但這超出了本模擬項目的范圍。通過以上五個部分的拆解你應(yīng)該已經(jīng)從“知道這個源碼包”變成了“能駕馭、能修改、能debug、能應(yīng)用到實際問題”的狀態(tài)。數(shù)字全息是一個深邃而有趣的領(lǐng)域這個像面離軸干涉的Matlab模擬是你探索計算成像世界一個非常扎實的起點。多動手改參數(shù)多觀察中間結(jié)果的變化你會對光與計算的融合有更直覺的理解。本文還有配套的精品資源點擊獲取