Extracting points from a matrix tangential to a centreline
3 次查看(过去 30 天)
显示 更早的评论
HB
2020-3-13
Hello,
I have a tube with two branches (please see attached point cloud showing this model where the columns are x y z and a variable respectively). I also have a centreline that runs through the main branch of the model. Please also see attached centreline points where the columns are x y z respectively. I have also attached images showing the model and centreline. The red outline drawn on the model geometry highlights where the main branch is located.
At each centreline point I am trying to extract the points located tangentially along the model to that centreline point to a new matrix along with their corresponding variable value. I have attached two images that hopefully clarify better what I want to achieve. I am pretty new to Matlab so in depth explanations would be great.
Thanks in advance.
6 个评论
darova
2020-3-13
The best idea i have: check each point of centerline (red) with data
ix = (xdata-xc(i)).^2 + (ydata-(yc(i)).^2 < R.^2; % indices inside circle
HB
2020-3-14
Thanks!
I tried the following script to check out your suggestion however I get the error "Array indices must be positive integers or logical values".
CL = load('Centerline.txt');
CLx = CL(:,1);
CLy = CL(:,2);
CLz = CL(:,3);
geom_variable=load('geom_variable.txt');
geomx = geom_variable(:,1);
geomy = geom_variable(:,2);
geomz = geom_variable(:,3);
ix = (geomx-CLx(i)).^2 + (geomy-(CLy(i)).^2 < R.^2;
HB
2020-3-15
Yes! Sorry. i = 46 and R= 1.5. ix matrix produced has the same number of rows as geomx, geomy and geomz. How do I then extract values that intersect with the radius ?
darova
2020-3-15
Maybe better be to use pdist2. IT calculates every possible combination of distances
You have two set of data (46 and 113502). So it creates matrix of size 46x11502
Try this
D = pdist2([CLx CLy CLz],[geomx geomy geomz]); % matrix of size 46x113502
D1 = D < 2; % distances smaller than "2" (matrix "0" and "1")
[~,ix] = find(sum(D1)>0); % if columns has one "1" - find it
plot3(CLx,CLy,CLz,'r') % centerline
hold on
plot3(geomx,geomy,geomz,'.b') % all data
plot3(geomx(ix),geomy(ix),geomz(ix),'.r') % data close to centerline
hold off
回答(2 个)
darova
2020-3-16
Here is what i did
I moved data to origin and rotated at some points of centerline
Once centerline is vertical i extract data -1<z<1 (for example)
I created surface at some points of centerline
countour was used to get crossection
Result
See attached script
3 个评论
HB
2020-3-17
Thanks, that is such a great way of doing it, very helpful!!
What do you think would be the best way of exporting the x,y,z values of the cross-sections to a text file and their corresponding variable value which are located in geom_variable(1:1:end,4)?
darova
2020-3-17
To write data
DATA = [];
for i = 1...
% code
V(:,1) = [];
DATA = [DATA V'];
end
dlmwrite('myFile.txt',DATA,'delimiter','\t')
- and their corresponding variable value which are located in geom_variable(1:1:end,4)?
Im afraid crossection line doesn't have corresponding values of geom_variable.
Closest points/values (red) can be found. But it's always different number
Change these lines in a loop
ix = abs(V(3,:)) < 0.1; % data for crossection creating
length(find(ix))
HB
2020-3-17
编辑:HB
2020-5-4
Ah I see, what I need is cross sections that consists of points from my geometry and the corresponding variable which is shear stress. The data is actually an artery model from a CFD.
What does the values length(find(ix)) values refer to? And how would I export the points shown in your latest plot?
Thanks
30 个评论
darova
2020-3-17
This line means the number of points
length(find(ix))
Once the number of points is different (very simple scheme)
cross1 cross2 cross3 ...
1 2 1
5 4 2
3 2
5
How do you want to have points? Separate txt?
darova
2020-3-18
this should work
ix = abs(V(3,:)) < 0.5; % data for crossection creating
dlmwrite([num2str(i) '.txt'],V(:,ix)','delimiter','\t')
HB
2020-3-18
It seems to have written one text file 45.txt containing one of the cross sections?
Also this cross section has 14 points, is there an easy way of extracting more points (as many as possible)?
Thanks
HB
2020-4-16
Great thanks, that works!! Is there any way to export the corresponding variable values which are in geom_variable(:,4)...
HB
2020-4-18
Thanks! Yes by corresponding variable, I meant that for every xyz position I have a variable value corresponding to that position.
Apologies for another question: I triedtthe script on another case (see attached files), but for some reason I am recieving the following error
"Error using linspace (line 22)
Inputs must be scalars.
Error in crossection (line 37)
tt = linspace(min(t)+0.02,max(t)-0.02,3000);"
Any suggestions why? Thanks!
HB
2020-4-28
Ah yes, you're correct, it is a unit issue!
Quick question - From the script I don't wholly understand how the the variable values are interpolated for each point? Do you mind clarifying? Thanks!
darova
2020-4-28
Imagine this is your data (x,y,z and c - some value - color)
Your data is like cylinder. I usually use interp2 or griddata for interpolation. Your data should like this (for each X Y only one Z)
So to make your data more 'flat' i unwided it (don't know if it's correct word). Look - >
[t0,r0] = cart2pol(x0,y0); % convert data from cartesian to polar
R = griddata(t0,z0,r0,T,Z); % calculate radius for grid points
C = griddata(t0,z0,c0,T,Z); % calculate color for grid points
[X,Y] = pol2cart(T,R); % convert back to cartesian
HB
2020-5-4
编辑:darova
2020-5-4
Thank you, that's a great explanation!
I was actually trying the script on a some other cases, and for some reason on a few of them points of the cross section line are missing as seen in image below. I don't know if there is something that can be modified in the script to fix this? I have attached the files FYI. Thanks!
HB
2020-5-4
Yes thanks, now there are no gaps! however for the variable value there are NaN values for some xyz positions..can these values be obtained/interpolated?
darova
2020-5-4
My bad. griddata has some problem at the boundaries. scatteredInterpolant uses extrapolation
I changed these lines
F = scatteredInterpolant(t(:),5*gz1(:),gc(ix));
gc1 = F(ct,ct*0);
HB
2020-5-7
Thanks, thats great!
I actually hope you don't mind me asking you for your opinion on something. I have second approach to achieving the above. Please see attached script. However it doesn't quite work for cases where there is a branch coming of the primary branch such as the cases I've shown you (because it also tries to extract the secondary branch). However it works fine when there is only the primary branch. It would be very helpful if you could maybe see if you can spot anything in the script that could be modified to make it work. If you don't think it can be, what are the issues with it in your opinion? Thanks a lot!
darova
2020-5-7
Sure, look on this picture
As i understood correctly: the script calculates red distance and if it's small enough the point is taken
You can calculate green distance too and compare it
d2 = pdist2(geomp(igeomp,1:3), CLp(iCLp,1:3)); % green distance
if dtmp<dist && d2 < 5 % compare red and green distance
Is it that you wanted?
HB
2020-5-7
编辑:HB
2020-5-7
Well the script extracts cross sections of the geometry at each CL point. So I think as well as the red distance it also calculated distances to the wall somehow.
If you were to run the attached case (no branches coming off primary) it does so fine. But not for the case I attached in the previous message, which I am hoping I can get to work with the script or figure out why it's not the best approach for cases with branches.
Your suggestions/assistance would be great! Thanks
darova
2020-5-7
- So I think as well as the red distance it also calculated distances to the wall somehow.
It doesn't that's why second branch is taken. I changed your script a bit
HB
2020-5-8
Thank you, that is so helpful!
I was wondering in your script how do you actually round off/close the gaps in the cross sections?
HB
2020-5-8
Ah sorry! I was referring to the script that you gave from my question on the 4th May. You edited the script to close gaps in the cross sections, I was just wondering how this was done? Also the script 'rounds off' the geometry where there is the branch to create the cross section, how is it done exactly? Is it the same way as the gaps were closed? Thanks!
darova
2020-5-8
HB
2020-7-8
编辑:darova
2020-7-8
Hello,
Hope you’re well!
I’ve been running cases with the latest script (see attached crossection_latest.m). It seems able to extract some cross sections but not others. For instance if I change line 26 from “i = 1:5:size(CL,1)-1” to “i = 1:size(CL,1)-1” therefore to include all the cross sections I get the error:
Error using linspace (line 22)
Inputs must be scalars.
Error in crossection_latest (line 43)
zz = linspace(min(gz1),max(gz1),10);
Any suggestions on how this can be overcome this (I have attached the files for a case)? It happens for a few of my cases, many thanks in advance!
darova
2020-7-8
Hello!
I think the problem you described is cause by radius you choosed. You didn't found all the points of an arteria
ix = sum(D<2) > 0; % if columns has one "1" - find it
You arteria is thicker than you think. Try to increase the radius
ix = sum(D<2.5) > 0; % if columns has one "1" - find it
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Live Scripts and Functions 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改,无法完成操作。请重新加载页面以查看其更新后的状态。
您也可以从以下列表中选择网站:
如何获得最佳网站性能
选择中国网站(中文或英文)以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)