Azzera filtri
Azzera filtri

decsg fails to decompose a few ellipses out of a set of 19 for unknown reason

19 visualizzazioni (ultimi 30 giorni)
Hello all, I am attempting to generate a decomposed geometry matrix of various ellipses in a square region using decsg in order to then extrude the geometry into 3D and mesh for use in a static structural model in the pde toolbox.
However, some of the ellipses (#4 when included with #s 1-3 or #19 when included with the rest) cause decsg to fail when included in the decsg inputs as so:
Error using decsg
Unable to decompose objects in the geometry description matrix into minimal regions.
Error in scatteredEllipses (line 42)
[dl, ~] = decsg(gd,sf,ns);
The scatteredEllipses function is presented below in its entirety, I have also attached the force_x, force_y, and force_z variables within forces.mat which are 512x512 double matrices of force magnitudes which represent desired ellipse locations.
clearvars -except force*; close all;
%% Fit ellipses to each traction island
tractionMask = abs(force_x) + abs(force_y) + abs(force_z); %combine x,y,z tractions to create mask of all nonzero elements
tractionMask(tractionMask ~= 0) = 1; %convert nonzero elements to 1
cc = bwconncomp(tractionMask); %get all connected ellipse regions
s = regionprops(cc,{'Centroid','Orientation','MajorAxisLength','MinorAxisLength'}); %fit ellipse parameters to each detected region
figure, imshow(tractionMask,[]), hold on
theta = linspace(0,2*pi);
orientationList = zeros(size(s,1),2);
for i = 1:length(s)
col = (s(i).MajorAxisLength/2)*cos(theta);
row = (s(i).MinorAxisLength/2)*sin(theta);
M = makehgtform('translate',[s(i).Centroid, 0],'zrotate',deg2rad(-1*s(i).Orientation));
[orientationList(i,1),orientationList(i,2)] = pol2cart(deg2rad(-1*s(i).Orientation),1);
D = M*[col;row;zeros(1,numel(row));ones(1,numel(row))];
plot(D(1,:),D(2,:),'r','LineWidth',2); hold on
s(i).Orientation = -1 * s(i).Orientation; %convert orientation to xy coordinate frame
%axis xy
hold off
%% Generate geometry
meshPtsFwdSol = 512;
halfSide = meshPtsFwdSol/2;
bound = [3; 4; -halfSide; halfSide; halfSide; -halfSide; halfSide; halfSide; -halfSide; -halfSide];
ns = char('bound');
sf = 'bound';
full = 1:length(s);
workingSubset = 4:length(s)-1;
g = 1;
for i = workingSubset %switch with full to cause error
sf = [sf,'+elli',num2str(g)];
ns = char(ns,['elli',num2str(g)]);
elli(:,g) = [4; s(i).Centroid(1)-halfSide; s(i).Centroid(2)-halfSide; s(i).MajorAxisLength/2; s(i).MinorAxisLength/2; deg2rad(s(i).Orientation); 0; 0; 0; 0];
g = g + 1;
gd = [bound,elli];
ns = ns';
[dl, ~] = decsg(gd,sf,ns);
pdem = createpde;
g = geometryFromEdges(pdem,dl);
facets = facetAnalyticGeometry(pdem,g,0);
gm = analyticToDiscrete(facets);
pdem.Geometry = gm;

Risposte (1)

Jaynik il 21 Mag 2024 alle 8:00
Modificato: Jaynik il 22 Mag 2024 alle 7:48
Hi Samuel,
After tinkering around with workingSubset, I found that ellipses #4 and #18 were causing the error when integrated with others and the code works for workingSubset = [1:3, 5;17, 19]. I started further analyzing #3 and #4 together by plotting the ellipses using the following code:
figure, hold on;
axis equal;
xlim([-halfSide, halfSide]);
ylim([-halfSide, halfSide]);
theta = linspace(0, 2*pi);
for idx = workingSubset
a = s(idx).MajorAxisLength/2; % Semi-major axis
b = s(idx).MinorAxisLength/2; % Semi-minor axis
X0 = s(idx).Centroid(1)-halfSide; % Adjusted X centroid
Y0 = s(idx).Centroid(2)-halfSide; % Adjusted Y centroid
phi = deg2rad(s(idx).Orientation); % Orientation in radians
% Parametric equation of the ellipse
x = X0 + a*cos(theta)*cos(phi) - b*sin(theta)*sin(phi);
y = Y0 + a*cos(theta)*sin(phi) + b*sin(theta)*cos(phi);
plot(x, y, 'r', 'LineWidth', 2);
hold off;
On plotting the ellipses, I could see that the ellipses #3 and #4 are very close to each other. This could be a cause for the error of invalid geometry.
To resolve this:
  • If possible, increase the distance between close geometries. If your geometries are very close to each other but not overlapping or intersecting, increasing the distance between them may help resolve the issue.
  • You can also try simplyfying the geometries. This could involve reducing the number of points in your geometries, or simplifying the shapes themselves.
According to the documentation of decsg, specifically this part:
Hope this help!




Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by