X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;ds=sidebyside;f=octave_packages%2Fgeometry-1.5.0%2Fgeom2d%2FcrackPattern.m;fp=octave_packages%2Fgeometry-1.5.0%2Fgeom2d%2FcrackPattern.m;h=3ec4ef5d1bbc217105abf37b6d6c28e446b73bb3;hb=c880e8788dfc484bf23ce13fa2787f2c6bca4863;hp=0000000000000000000000000000000000000000;hpb=1705066eceaaea976f010f669ce8e972f3734b05;p=CreaPhase.git diff --git a/octave_packages/geometry-1.5.0/geom2d/crackPattern.m b/octave_packages/geometry-1.5.0/geom2d/crackPattern.m new file mode 100644 index 0000000..3ec4ef5 --- /dev/null +++ b/octave_packages/geometry-1.5.0/geom2d/crackPattern.m @@ -0,0 +1,189 @@ +%% Copyright (c) 2011, INRA +%% 2007-2011, David Legland +%% 2011 Adapted to Octave by Juan Pablo Carbajal +%% +%% All rights reserved. +%% (simplified BSD License) +%% +%% Redistribution and use in source and binary forms, with or without +%% modification, are permitted provided that the following conditions are met: +%% +%% 1. Redistributions of source code must retain the above copyright notice, this +%% list of conditions and the following disclaimer. +%% +%% 2. Redistributions in binary form must reproduce the above copyright notice, +%% this list of conditions and the following disclaimer in the documentation +%% and/or other materials provided with the distribution. +%% +%% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +%% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +%% IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +%% ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE +%% LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +%% CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +%% SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +%% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +%% CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +%% ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +%% POSSIBILITY OF SUCH DAMAGE. +%% +%% The views and conclusions contained in the software and documentation are +%% those of the authors and should not be interpreted as representing official +%% policies, either expressed or implied, of copyright holder. + +%% -*- texinfo -*- +%% @deftypefn {Function File} {@var{e} = } crackPattern (@var{box}, @var{points}, @var{alpha}) +%% Create a (bounded) crack pattern tessellation +%% +%% E = crackPattern2(BOX, POINTS, ALPHA) +%% create a crack propagation pattern wit following parameters : +%% - pattern is bounded by area BOX which is a polygon. +%% - each crack originates from points given in POINTS +%% - directions of each crack is given by a [NxM] array ALPHA, where M is +%% the number of rays emanating from each seed/ +%% - a crack stop when it reaches another already created crack. +%% - all cracks stop when they reach the border of the frame, given by box +%% (a serie of 4 points). +%% The result is a collection of edges, in the form [x1 y1 x2 y2]. +%% +%% E = crackPattern2(BOX, POINTS, ALPHA, SPEED) +%% Also specify speed of propagation of each crack. +%% +%% +%% See the result with : +%% figure; +%% drawEdge(E); +%% +%% @seealso{drawEdge} +%% @end deftypefn + +function edges = crackPattern(box, points, alpha, varargin) + + if ~isempty(varargin) + speed = varargin{1}; + else + speed = ones(size(points, 1), 1); + end + + % Compute line equations for each initial crack. + % The two 'Inf' at the end correspond to the position of the limit. + % If an intersection point is found with another line, but whose position + % is after this value, this means that another crack stopped it before it + % reach the intersection point. + % There is one 'end position' for each side of the crack. + lines = [points speed.*cos(alpha) speed.*sin(alpha) Inf*ones(size(points, 1), 2)]; + + % initialize lines for borders, but assign a very high speed, to be sure + % borders will stop all cracks. + dx = (box([2 3 4 1],1)-box([1 2 3 4],1))*max(speed)*5; + dy = (box([2 3 4 1],2)-box([1 2 3 4],2))*max(speed)*5; + + % add borders to the lines set + lines = [lines ; createLine(box, dx, dy) Inf*ones(4,2)]; + + edges = zeros(0, 4); + + + while true + modif = 0; + + % try to update each line + for i=1:size(points, 1) + + % compute intersections with all other lines + pi = intersectLines(lines(i,:), lines); + + % compute position of all intersection points on the current line + pos = linePosition(pi, lines(i,:)); + + % consider points to the right (positive position), and sort them + indr = find(pos>=0 & pos~=Inf); + [posr, indr2] = sort(pos(indr)); + + + % look for the closest intersection to the right + for i2=1:length(indr2) + + % index of intersected line + il = indr(indr2(i2)); + + % position of point relative to intersected line + pos2 = linePosition(pi(il, :), lines(il, :)); + + % depending on the sign of position, tests if the line2 can + % stop the current line, or if it was stopped before + if pos2>0 + if pos20 + if pos2