發表文章

(幾何)相機標定_投影轉換_針孔成像_相機內外部參數_齊次座標(homogeneous)

圖片
  (幾何) 相機標定(校正)/( Geometric ) Camera calibration 是用來求初將三維世界投影成二維影像所需要的各種參數,這些參數也就是攝影機參數。 也可理解為從世界座標系換到二維影像座標系的過程,也就是求出最終投影矩陣的過程,用於估計成像過程的數學參數。 摘自 Medium-Camera Calibration 相機校正 小(針)孔成像(Pinhole) 在進入到相機標定前,需要先掌握針孔成像的基本原理。 以下是一張圖,從右側去看是一個世界坐標系中的一棟房子,左側是一個成像的螢幕。 如何在成像螢幕上去呈現出清晰完整的房子影像呢? 摘自youtube影片截圖- Pinhole and Perspective Projection | Image Formation 最簡單方法是透過針孔,開一個小孔的薄板。 摘自 https://kknews.cc/zh-mo/other/mlkmga2.html 以前國中時候學過的自然科學實驗,會拿蠟燭發光的光線通過針孔後,投影在後方紙屏上,蠟燭燭光在照射到後面紙屏之前,在中間額外再放一個隔板一個小孔,可觀察到上下顛倒左右相反的實像。 如果目標成像紙屏往後移動遠離中間隔板,則會發現針孔成像會放大,影像相對變比較暗。 如果目標成像紙屏往前移動靠近中間隔板,則會發現針孔成像會縮小,影像相對變比較亮。 此外還滿足以下比例關係 \[ \frac{\text{物長}}{\text{像長}} = \frac{\text{物距}}{\text{像距}} \] 回過頭去看Perspective Imaging with Pinhole 的圖繼續去探討,中間虛線是光軸也就是垂直於影像平面的軸,在空間座標戲中建構三維xyz軸座標。而針孔到影像平面的距離就是有效焦距 f。 真實的三維空間以房子頂部的一個點 \(P_o\),它在「相機座標系」中的三維座標寫成: \[ \mathbf r_o=(x_o,y_o,z_o) \] 其中 \(x_o\) 是左右位置、\(y_o\) 是上下位置、\(z_o\) 是物體離相機沿光軸方向的深度。 光線從物體點 \(P_o\) 出發,經過 Pinhole,最後與左邊的成像平面相交,得到影像點 \(P_i\)。 因此: \[ P_o \rightarrow \text{Pinhole} \...

Head Pose estimation_PnP(Perspective-n-Point)多點透視問題_ 2D and 3D feature-based alignment

圖片
頭部姿態的三個自由度 A survey of head pose estimation methods 在2020年一篇IEEE論文是有關於「頭部姿態估計方法綜述」的回顧中,可得知頭部姿態通常以三個自由度來表示。分別是俯仰角(pitch)、偏航角(yaw)與翻滾角(roll),頭部姿態也蘊藏豐富的資訊,比方點頭可能代表答應、理解等,搖頭則表示否定。 延伸應用時常用於視線估計、臉部表情分析,於駕駛監控系統中也十分重要,能去評估駕駛精神狀況。也有應用是分析學生上課專心程度。 傳統頭部姿態估測方法可先偵測人臉特徵點,再藉由頭部模型建立二維影像特徵點與三維模型點之間的對應關係,以估測三維頭部姿態。  一些先備知識可以再去這篇溫故 (幾何)相機標定_投影轉換_針孔成像_相機內外部參數_齊次座標 特徵式對齊 feature-based alignment (特徵式對齊)是估計兩組或多組已匹配的 2D 或 3D points之間運動關係的問題。 在feature-based alignment中,有一種相當常見的特定情況,就是根據一組 2D 點的投影位置,估計物體的 3D 姿態。關於姿態估計問題也被稱為是所謂的外部參數校正。與其相對的是相機內部的參數校正,比方焦距參數。 從三個對應點中恢復姿態,需要的信息是最少的,稱為「透視三點問題」(perspective-3-point-problem, P3P)。當擴展到更多點的問題,合起來稱為「PnP」。 直接線性變換 (direct linear transform, DLT) 摘自 Computer Vision Algorithms and Applications, Richard Szeliski, 2010 - Chapter 6 Feature-based alignment-6.2 Pose estimation 已知p(u,v)、p(X,Y,Z),求出K、R、t。 我們需要知道一組2D-3D對應,並知道對應點的2D與3D座標。 將上面式子中的K[R|t]展開,變成三列四行的矩陣,K內參矩陣由於固定可省略不影響推導。 兩側相乘可得到如下三個方程 其中最底下的z 會等於右側的最底下的式子 我們可把小z帶入上面式子的zu 、 zv去 化簡變底下式子 我們的未知就是L 11到 L 34這些參...

DSP筆記6_功率頻譜密度 Power Spectral Density(PSD)_Periodogram(週期圖)_用scipy內建的periodogram跟最原始簡單的作法比較

圖片
事實上,在上一篇我們完全不依賴scipy手刻的PSD計算的python程式 它的演算法流程就是所謂的Periodogram(週期圖)最原始簡單的作法。 periodogram(週期圖)是一種用來估計訊號功率頻譜密度 使用單一區段全部資料來計算 PSD 雖然速度快,但具有對雜訊敏感且變異數較大的限制。 顯示訊號的功率如何分布在不同頻率上。 常用於signal processing、vibration analysis、PPG/EEG/ECG 分析 在數學上,對於一個長度為  N  的離散時間訊號  x [ n ] : 會用以下數學式子表示週期圖功率頻譜密度 推導過程 Step1.計算離散傅立葉轉換(實際應用場景中,可用快速傅立葉轉換更有效率做出來。) Step2.計算功率譜,也就是計算「振幅(幅度大小)平方」 Step3.除以樣本數N做正規化 原始periodogram功率頻譜密度 提供 每個樣本的功率估計值 。 原始 periodogram就是 最初版本、最簡單的 PSD 估計方法。 在實務上,通常使用 FFT 來有效率地進行計算。 原始periodogram功率頻譜密度特性 變異數 : 高;即使增加 \(N\),變異數也不會下降,因此結果通常會比較多雜訊干擾。 偏差(Bias)、頻譜洩漏(Spectral Leakage): 在有限長度訊號中的不連續現象,會使功率擴散到鄰近的頻率。 解析度(Resolution) Δ f = f s N ;較長的訊號可以提升解析度。 一致性(Consistency): 原始 periodogram是不一致的估計量 原始periodogram功率頻譜密度限制 High variance : 峰值會產生波動 Spectral leakage : 有限長度資料會造成頻譜能量向周圍頻率擴散。 Poor consistency : N→∞時,原始 periodogram 不會收斂。 Trade-off between bias and resolution : 加窗(Windowing)可以降低頻譜洩漏,但會使頻譜峰值變寬。 常見改善方法 加窗(Windowing): Hann、Hamming、Blackman 視窗可降低 頻譜洩漏(spectral leakage) 。 Welch 方法(Welch'...

DSP筆記5_功率頻譜密度 Power Spectral Density(PSD)_從時域訊號到幅度大小的平方再除以資料總長度

圖片
功率頻譜密度(Power Spectral Density, PSD) 用來描述訊號在不同頻率下的功率分布情形。 功率頻譜密度(Power Spectral Density, PSD) S x ( f ) S_x(f) 適用於 功率訊號(power signals) 以及 廣義平穩程序(wide-sense stationary processes) 。 定義, 透過自相關函數(Autocorrelation)根據( Wiener–Khinchin theorem,維納–辛欽定理 )可得到如下式子 另一種定義 (此公式十分重要) : 總功率與性質(Total Power and Properties) 換言之,總功率(Total Power)可以由功率頻譜密度 S x ( f ) S_x(f) 對所有頻率積分後得到。 把給定區間的功率譜密度全部加總起來 練習程式ver1.建立一個典型的10Hz正弦波訊號 import numpy as np import matplotlib.pyplot as plt # Sinusoid parameters fs = 100 # Sampling freq (Hz)每秒取樣 100 次。 t = np.arange( 0 , 1 , 1 /fs)#[0. 0.01 0.02 ... 0.98 0.99],每次增加 0.01。 f0 = 10 # signal frequency # Sinusoidal signal x = np.sin( 2 * np.pi * f0 * t) # Plot plt.figure(figsize=( 8 , 4 )) plt.plot(t, x) plt.title( "10 Hz Sinusoid" ) plt.xlabel( "Time (s)" ) plt.ylabel( "Amplitude" ) plt.grid( True ) plt.show() 現在產生的是: x ( t ) = sin ⁡ ( 2 π f 0 t ) x(t)=\sin(2\pi f_0t) 其中 f 0 = 10  Hz f_0=10\text{ Hz} 意思是一秒鐘會完成...

DSP筆記4_能量頻譜密度 Energy Spectral Density(ESD)_從時域訊號到幅度大小的平方

圖片
  能量頻譜密度 Energy Spectral Density(ESD) 用於描述訊號在不同頻率下能量分布情形。 用  S E ( f )  來表示,X(f)為原先時間域x(t)的訊號做傅立葉轉換的結果。 能量譜密度的總能量關係滿足如下式子,就稱作帕塞瓦爾定里Parseval's theorem)。 總能量會等同於對時域跟頻域函數從負無窮到正無窮取平方的積分 所以 練習程式ver1.矩形脈衝Energy Signal import numpy as np import matplotlib.pyplot as plt fs = 1000 #Sampling Frequency T = 1 # Duration (seconds) t = np.linspace( 0 , T, int (T*fs), endpoint= False ) #建立一個全為0的訊號 x = np.zeros( len (t)) x[ int ( 0.4 *fs): int ( 0.6 *fs)] = 1 # 0.4 ~ 0.6 秒設定為 1 # Plot plt.figure(figsize=( 8 , 4 )) plt.plot(t, x) plt.title( "Rectangular Pulse" ) plt.xlabel( "Time (s)" ) plt.ylabel( "Amplitude" ) plt.grid( True ) plt.show() 運行效果 可以先得知目前此時間訊號的函數區間定義如下 x ( t ) = { 1 , 0.4 ≤ t < 0.6 0 , 其他時間 x(t)= \begin{cases} 1, & 0.4\leq t<0.6\\ 0, & \text{其他時間} \end{cases} 也就是一個寬度 0.6 − 0.4 = 0.2  秒 0.6-0.4=0.2\text{ 秒} 的矩形脈衝,是一個典型的 Energy Signal。。 它的能量為: E = ∫ − ∞ ∞ ∣ x ( t ) ∣ 2 d t E=\int_{-\infty}^{\infty}|x(t)|^2dt ...