not a number

安安
  • About
  • Publications
  • Archive
  • feeds

Posts match “ Matlab ” tag:

almost 9 years ago

將 Matlab 矩陣匯出成 Amira mesh (.am)

簡單筆記一下今天工作的成果:將 Matlab 的矩陣噴出成 Amira 的 .am 檔案,也就是 Amira 3D mesh。
這個成果是根據我手上有的 .am 檔案還有 Matlab Central 上的 LoadData_Amira.m 逆推的,特此致謝。

雖然從 LoadData_Amira.m 看起來 AM 檔案有支援 uint8 之外的其他數值格式,但是我剛才弄半天弄不太出來,所以為了簡便,就先強制輸出為 uint8 吧,科科!

Read on →
  • Matlab
  • .m
  • Amira
  • confocal
  • August 30, 2013 18:01
  • Permalink
  • Comments
 
almost 9 years ago

把 Matlab 裡面 2D 或者更高維度的矩陣當作 vector 噴出來

剛才突然想到,之前一直有個很差的習慣,就是遇到2D的矩陣,就用很彆扭的方式去解決,例如要解2D矩陣的 max ,就寫

A = rand(6, 3);
max_num = max(max(A));

照理來說min和mean這類 會從 column 開始做起然後才做更高維度 的函數,我也都是用一樣的方法來解決。

可是俗話說的好, 夜路走多了會遇到鬼,如果遇到 std 我就沒轍了,如果你這樣寫你就準備出糗了

Read on →
  • Matlab
  • matrix
  • Array
  • vector
  • 3D
  • 2D
  • September 10, 2013 14:09
  • Permalink
  • Comments
 
almost 9 years ago

自訂 Matlab polar 函數的邊界

有在用 Matlab 的人可能有用過 polar這個畫極座標的繪圖函數,但是很令人崩擴的是 polar 這個函數本身功能非常有限,如果要多別的事情,都要靠 handler 還有 set搭配來完成,不像 plot早就已經做到萬能的地步,幾乎想要做什麼,都可以直接在 argument 裡面完成。

現在的問題是,假設我想要把兩個函數都畫在同一張 polar 上面,聰明的各位可能馬上就會想到可以透過 hold on,來實作,這是沒錯的,但是如果實際去做,會發現一個很糟糕的事情:polar的limit會被第一個畫的函數所限制住,例如假設我想要畫以下的函數:

最簡單的寫法大概是這樣:

close all;

N = 101;
theta = linspace(0, 2*pi, N);

y1 = abs(cos(theta));
y2 = 1.5 * abs(sin(theta));

figure(1),
polar(theta, y1, 'k');
hold on,
polar(theta, y2, 'r');
hold off;

但是畫出來卻會變這樣:

囧!因為 limit 已經被第一個polar決定了,所以如果後面再hold on上去的圖案更大,就只能悲劇。

Read on →
  • Matlab
  • plot
  • polar
  • 論文
  • September 26, 2013 10:58
  • Permalink
  • Comments
 
over 8 years ago

解決Amira不能開Matlab的mat檔

Amira 是一個許多做顯微術或者生物影像的人都會用到的軟體,可以用它來具象化許多3D的data。
除了原生的AM檔案之外,有時候我們也會希望藉由 Amira來顯現來自 Matlab我們加工後或者計算後的檔案,因此就會需要由 Amira來打開 .mat檔案;照理來說 Amira 是支援的,但是我自己實際測試,卻會噴出這樣的錯誤訊息:

Reading minivirus_A.mat ...
opening connection to Matlab engine...
Couldn't open libeng.dll:
File does not exist.
Couldn't open libmat.dll:
File does not exist.
Couldn't open libmx.dll:
File does not exist.
Failed to connect to Matlab
Be sure Matlab is correctly installed, and the path is valid
Error opening file G:/minivirus_A.mat

網路上查半天,都找不到問題所在或者解法,但是自己嘗試去解決,卻發現解法超乎想像地簡單。

首先先試著去找 libeng.dll 、libmat.dll或libmx.dll這三個鬼東西在哪。答案很簡單,就在

C:\Program Files\MATLAB\R2012b\bin\win64

的下面。我想 Amira 抓不到的原因非常簡單,是因為 Amira 不知道 64-bit 的 Matlab 已經改變存放這些lib的位置了。

解決方法

解法非常簡單,就是改系統變數中的 PATH值。

首先對著 我的電腦 按右鍵,選 內容 ,接著選 進階系統設定 ,然後按 進階 這個頁籤,最後按最下面的按鈕 環境變數 。接著在下面的小框框中找到PATH,點選 編輯,然後在最後面加上

;C:\Program Files\MATLAB\R2012b\bin\win64

就大功告成,此時重開 Amira 後就可以開啟 .mat 檔案了。

  • Matlab
  • Amira
  • path
  • November 28, 2013 13:41
  • Permalink
  • Comments
 
about 8 years ago

求圓球的質量投影

之前就有想過要用解析解求一個均勻的圓球的質量投影,但是推半天沒推出來,只好先用MATLAB做出一個3D的volume,再用 sum(V,3)來求投影質量。我知道這樣超白癡的,而且要是球稍微大一點,電腦就記憶體不足當掉給你看了。

剛才突然... 就夢到正確的解析解做法了......

原理

已知圓球球面的公式是\

其中 為一常數,是球面的半徑。

現在我們將其推廣,把球面推廣成球體,則參數 要滿足

其中

所以要快速求出一個球體的質量投影,只要移項還有乘以2就好了,得到:

不過記得要把超出 support 的地方手動設為零,不然可是會出現 complex 的。請直接參考下面的實作程式碼。

MATLAB 實作

sphere_proj.m
function Z = sphere_proj(DIA)

RAD = DIA/2; % radius

x = linspace(-RAD, RAD, DIA);
[X, Y] = meshgrid(x, x);
R = sqrt( X.^2 + Y.^2);
Z = 2 * sqrt( RAD^2 - R.^2 );
Z( R > RAD) = 0;
  • Matlab
  • June 17, 2014 14:40
  • Permalink
  • Comments
 
over 7 years ago

EEMD (Ensemble Empirical Mode Decomposition) 訓練講義

簡介

訓練目標

以1D為切入點,了解 EEMD 之基本原理,並以 Matlab 語言實作之。

Matlab

完成這份訓練講義之後,預期學習者將會應用到或者學到如何

  • 創造並使用 function
  • 找出資料向量中的極值
  • 使用 interp1 內插出平滑的資料曲線
  • 打造 EMD 引擎
  • 將 white noise 加入資料之中
  • 完成 EEMD 函數

習題

尋找極值函數 findExtrema

在這題當中,我們要學會寫出一個傳入資料向量 data,可以回傳資料中極值的 值 以及 位置 的函數。以下將分為四個小題進行。

學會如何自定函數

試找尋資料,學習如何創造一個自定函數,或查詢 Matlab 內建教學

doc function

或參考以下範例:

numAdd.m
function y = numAdd(a, b)
y = a + b;

尋找最大值

試寫一個函數 findMax.m,輸入任意長度向量 data 可傳回兩組輸出:

  • 最大值所在的位置 max_ind
  • 最大值的值 max_val

提示: 函數的第一行可能長成這個樣子:

findMax.m
function [max_ind, max_val] = findMax(data)
% coding here

請注意,一串資料可以有多個 local maxima。此題要求的是找出所有的 local maxima,而非僅是資料中的最大值,否則用內建函數 max() 即可達成。

尋找最小值

承上題,找出「所有」最小值 min_val 及其分別之所在位置 min_ind。

提示: 試想,是不是將上一題的傳入資料由 data 改為 -data 就可以達成目標了?

尋找所有極值

綜合前面三題所學,請完成一個可以找出所有極值,並分別回傳最大最小值的值與位置的函數 findExtrema.m

findExtrema.m
function [max_ind, max_val, min_ind, min_val] = findExtrema(data)
% coding here

注意: data(1) 與 data(end) 是否也須考慮為「極值」?如果是,其原因為何?

內插曲線

下面這兩題,我們將會利用前一大題所得到的資料,以及 Matlab 內建的 interp1函數,得到 EMD 所需要的「趨勢」曲線。

使用離散資料內插曲線

利用 interp1,將 max_ind 和 max_val內插成一條長度和取樣率皆與 data 相同之「連續」曲線 max_trend。
利用 interp1,將 min_ind 和 min_val內插成一條長度和取樣率皆與 data 相同之「連續」曲線 min_trend。

注意: max_trend或min_trend的頭尾,都要和 end 一樣(可以想成是 boundary conditions),否則會有些區域沒有定義。

求得「平均」曲線

將上小題的兩條曲線做點對點平均,得到 avg_trend。

打造 EMD 引擎

參考 paper 上的演算法,以及上面所練習的技巧與成品,應該已經可以完成 EMD 的引擎。其界面如下:

emd.m
function modes = emd(data, nmode)
% coding here

其中 nmode 是所要求的模態數量。請注意,如果原本的 data 長度是 的話,則得到的輸出資料 modes的大小會是 。

打造 EEMD 引擎

如果前面的 emd.m 已經驗證完成,那麼接下來就是讓它進化成 EEMD 的時候了。比起 EMD , EEMD 多了下面兩項特色

  • 外加了額外的 white noise 進去
  • 先跑出 N 組各別使用不同的 white noise 所解出來的結果,再一起取平均。

下面我們只拆成兩小題練習,不足之處請直接參閱 paper。

加入 white noise

假定使用者指定一個特定的 white noise 強度 wn_std,我們有兩種方式可以將等效強度的 white noise Wn 加到 data 裡面去,第一種是:

其中 指 data, 是求標準差的函數;或者第二種方法

試問:兩種方法何種較優?較劣者,有什麼致命的缺點?

請利用 Matlab 內建的 randn() 完成此題。

完成 EEMD 架構

參考 paper 上的演算法,以及上面所練習的技巧與成品,終於我們可以完成 EEMD 的全面架構。完成品的輸入資料除了 data以及nmode之外,還多了指定 ensemble 數量的 nens 與 white noise 的標準化強度 wn_std。其界面如下:

eemd.m
function modes = eemd(data, nmode, nens, wn_std)
% coding here
  • Matlab
  • EMD
  • EEMD
  • October 14, 2014 18:52
  • Permalink
  • Comments
 

Copyright © 2013 leeneil . Powered by Logdown.
Based on work at subtlepatterns.com.