- Apply the technique to the inner region of your volume, far enough away from the edges that you won't have problems. Then the output matrix will be smaller than the original input.
- Write code to handle the edges as a special case. This can get complicated, and it is not obvious how exactly you should handle the special case. If you don't compute the contributions to the convolution that come from off-edge, it is like assuming the off-edge values are zero, which can give bad results for derivatives.
- Extend the matrix beyond the edges with an "upside down and backwards" extension. See 1-D example below. Then you filter the extended matrix. Then you keep only the central region of the filtered extended signal. This is what I would do. This is used for derivatives and other filtering of 1-D signals.
A non centered method to compute velocity derivative
5 ビュー (過去 30 日間)
古いコメントを表示
I have three vector vx, vy and vz, these vectors contain the values of velocity in every point of my 3d domain and in every time frame. I have a 3d domain and I'm computing derivatives of the velocity in every direction with the Sobel Operators in this way:
s1(:,:,1) = [1 2 1; 0 0 0; -1 -2 -1];
s1(:,:,2) = [2 4 2; 0 0 0; -2 -4 -2];
s1(:,:,3) = [1 2 1; 0 0 0; -1 -2 -1];
s2(:,:,1) = [1 0 -1; 2 0 -2; 1 0 -1];
s2(:,:,2) = [2 0 -2; 4 0 -4; 2 0 -2];
s2(:,:,3) = [1 0 -1; 2 0 -2; 1 0 -1];
s3(:,:,1) = [1 2 1; 2 4 2; 1 2 1];
s3(:,:,3) = -[1 2 1; 2 4 2; 1 2 1];
s1 = s1./32;
s2 = s2./32;
s3 = s3./32;
sx = s2;
sy = s1;
sz = s3;
dvxdx = imfilter(vx,sx,'conv'); Derivative of the x-component of the velocity along direction x
dvxdy = imfilter(vx,sy,'conv'); Derivative of the x-component of the velocity along direction y
dvxdz = imfilter(vx,sz,'conv'); Derivative of the x-component of the velocity along direction z
Since this method is a centered method, it doesn't work very well when I want to compute the derivatives on the boundary of my domain. Does anyone know a not centered method to compute derivatives of the velocity along all the directions given the three vectors vx, vy, vz?
0 件のコメント
採用された回答
William Rose
2022 年 1 月 19 日
編集済み: William Rose
2022 年 1 月 20 日
[edited to crrect typographical errors]
If you use a non-centered-to-the-right method, you will have issues on the right edge. If you use non-centered-to-the-left method, you will have issue on the left edge. You have several options:
Example: Extend a 1D signal by m on each side, so that it can be filtered with a (2m+1)-point-wide centered filter:
t=(0:20)/20; x=cos(2*pi*t)-2*t; %any waveform or recording
m=3; %number of points to add on each end
n=length(x);
%next: compute initial value of extended x
xext=[zeros(1,m),x,zeros(1,m)];
%next: compute the upside-down-and-backwards extensions
for i=1:m
xext(i)=2*x(1)-x(m+2-i); %before the beginning
xext(n+m+i)=2*x(end)-x(end-i); %tail: "put your behind in your past"-Puumba
end
text=(-m:20+m)/20; %extended time vector, for plotting
plot(t,x,'-b.',text,xext,'-rx'); legend('original','extended'); grid on
You can extend this idea to 3D.
5 件のコメント
William Rose
2022 年 1 月 26 日
I'd call it the Rose method! The face extension method is a generalizaiton of what has been done for 1D signals, in various papers. The methods for the edges and corners are original. You're welcome. Good luck with your work.
その他の回答 (0 件)
参考
カテゴリ
Help Center および File Exchange で Surface and Mesh Plots についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!