calcarea.m

这个GUI程序让你测量面积和图像中路径的长度

% hObject    handle to run (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hint: get(hObject,'Value') returns toggle state of run



% hObject    handle to calc (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

    function valid=testIfCanCloseATriangle(n,m,k,connectionsMat,points,triangleCenters)
        % tests if the three close a rectangle inside the shape,
        % if the new line doesn't cross any existing line
        % if the new triangle doesn't contain in it any existing triangle
        valid = false;
        
        angle1 = mod( angle(( points(k,:)-points(n,:) )*[1 ;i])+2*pi, 2*pi);
        angle2 = mod( angle(( points(k,:)-points(m,:) )*[1 ;i])+2*pi, 2*pi);
        middleAngle = (angle1+angle2)/2;
        if max([angle1,angle2] - middleAngle) > pi/2
            % The middle angle croses the wide angle... Neads to flip it.
            middleAngle = mod(middleAngle+pi,2*pi);
        end
        if min( abs(middleAngle - [angle1,angle2]-pi/2) ) < 1e-1
            % The three points are on a straight line
            valid = false;
            return
        end
        dAngle = min( abs([angle2 angle1] - middleAngle) );
        p1 = points(k,:);
        
        % Deciding if the rectangle is in the overall shape or not
        % The middle Angle line should cross odd number of exterior lines
        [a b] = find(connectionsMat == 1);
        f = find(b > a & b~=k & a ~=k);    % working only on half the matrix (not to count same line twice)
        a = a(f);   b=b(f);
        decision = false;
        itt = 0;
        nTests = 5;
        crossings = zeros(length(f),nTests);
        while ~decision && itt < 10
            for mm =1 :nTests 
                % testing for nTests number of lines coming out of the head
                % point. Comes to verify I don't acediently cross exactly on a
                % junction point.
                testAngle = middleAngle - dAngle + 2*dAngle*rand(1,1);
                p2 = points(k,:)+[cos(testAngle) sin(testAngle)]*-1e4; %don't know why the minus
                for nn=1 : length(f)
                    crossings(nn,mm) = doesLinesCross(p1,p2, points(a(nn),:), points(b(nn),:),'' );
                end        
            end
            itt = itt+1;
            crossingCount = sum( crossings );
            isPointIn = mod(crossingCount,2);

            %% decision is only true if all lines returned the same result (even
            %% or odd number of crossings)
            decision = true;
            if  sum(abs(diff(isPointIn)))
                decision = false;
            end
        end %while ~decision & itt < 10
        %%%%%%%%%%%%%%%%

        if ~( isPointIn(1,1) )
            % There where even number of crossings
            valid =false;
            return
        end
        
        % Veryfing that the new line doesn't cross any existing line
        % NEADS TO VERIFY THIS !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
        cross = lineCrossExistingLines(m,n,connectionsMat,points);
        if cross
            valid = false;
            return
        end
        
        %% Verifying the triangle doesn't contain any other triangle in it
        cross = zeros(3,1);
        for nn=1: size(triangleCenters,1)
            p1 = triangleCenters(nn,:);
            p2 = rand(1,2).*[1e4 1e4];
            cross(1) = doesLinesCross( p1,p2, points(m,:), points(n,:),'' );
            cross(2) = doesLinesCross( p1,p2, points(m,:), points(k,:), '' );
            cross(3) = doesLinesCross( p1,p2, points(k,:), points(n,:), '' );
            if mod( sum( cross),2 )
                valid = false;
                return;
            end
        end
        %%%%%
        valid = true;
                
        
        
    function cross = lineCrossExistingLines(p1,p2,connectionsMat,points)
        [a b] = find(connectionsMat);
        f = find( b > a &  ~(b==p1 & a==p2) & ~(a==p1 & b==p2) );
        N = length(f);
        a = a(f);   b=b(f);
        nn=0;
        cross = 0;
        while nn < N
            nn = nn+1;            
            cross = doesLinesCross( points(p1,:), points(p2,:), points(a(nn),:), points(b(nn),:),'EdgeCross' );
            if cross
                return
            end
        end
        
        
    function cross = doesLinesCross(pA1,pA2,pB1,pB2,mode)  
%         cross = 0     - no crossing
%         cross = 1     - there is a crossing
        cross = 0;
        a1Sat = false;  a2Sat = false;
        d_a = pA2-pA1;
        d_b = pB2-pB1;
        a1 = tan( angle( d_a(1)+i*d_a(2) ) );
        if abs(a1) > 1e5
            a1Sat = true;
        end
        b1 = pA1(2)-a1*pA1(1);
        a2 = tan( angle( d_b(1)+i*d_b(2) ) );
        b2 = pB1(2)-a2*pB1(1);
        if abs(a2) > 1e5
            a2Sat = true;
        end
        if a2Sat
            if a1Sat
                % both saturated
                if ~ ( min(pA1(2),pA2(2)) >= max(pB1(2),pB2(2)) || max(pA1(2),pA2(2)) <= min(pB1(2),pB2(2)) ) && abs( pA1(1)-pB1(1) ) < 0.1
                    cross = 1;
                    return
                end
            else
                % only a2 is sat
                x0 = (pB1(1)+pB2(1))/2;
                y0 = a1*x0+b1;
                if y0 >  min(pA1(2),pA2(2)) && y0 < max(pA1(2),pA2(2)) && ...
                    x0 > min(pA1(1),pA2(1)) && x0 < max(pA1(1),pA2(1)) && ...
                    y0 > min(pB1(2),pB2(2)) && y0 < max(pB1(2),pB2(2))
                    cross = 1;
                end
            end
        elseif a1Sat
                % only a1 is sat
                x0 = (pA1(1)+pA2(1))/2;
                y0 = a2*x0+b2;
                if y0 >  min(pA1(2),pA2(2)) && y0 < max(pA1(2),pA2(2)) && ...
                    y0 > min(pB1(2),pB2(2)) && y0 < max(pB1(2),pB2(2)) &&...
                     x0 > min(pB1(1),pB2(1)) && x0 < max(pB1(1),pB2(1))
                    cross = 1;
                end
        else
            if a1 ~= a2
                x0 = (b2-b1)/(a1-a2);
                y0 = a1*x0+b1;
                if x0 < max(pA2(1),pA1(1))-eps && x0 > min(pA2(1),pA1(1))+eps && x0 < max(pB2(1),pB1(1))-eps && x0 > min(pB2(1),pB1(1))+eps
                     cross =1;
                end
            elseif ~(max(pA1(1),pA2(1)) < min(pB1(1),pB2(1)) || min(pA1(1),pA2(1)) > max(pB1(1),pB2(1)) )
                cross =1;
                return
             end
        end
        
        if strcmp(mode,'EdgeCross')  && cross
            x = [pA1(1), pA2(1), pB1(1), pB2(1)];
            y = [pA1(2), pA2(2), pB1(2), pB2(2)];
            if min( (x-x0).^2 + (y-y0).^2 ) < 0.001
                cross =0;
            end
        end
        
function area = calcAreaMonteCarlo(points,connectionsMat)    
    % monte carlo !
    
    Y = [ min(points(:,2)) max(points(:,2)) ];
    X = [ min(points(:,1)) max(points(:,1)) ];
    dx = X(2)-X(1);   dy = Y(2)-Y(1);
    
    count = 0;
    N = 1e4;
    n =0;
    hWaitBar = waitbar(0,'Monte Carloing','CreateCancelBtn','closereq ');
    while n < N && ishandle(hWaitBar);
        p1 = rand(1,2) .* [dx dy] + [X(1) Y(1)];
        p2 = rand(1,2).*[1e4 1e4];
        
        [a b] = find(connectionsMat == 1);
        f = find(b > a);    % working only on half the matrix (not to count same line twice)
        a = a(f);   b=b(f);
        crossings = zeros(length(f),1);
        for nn=1 : length(f)
            crossings(nn) = doesLinesCross(p1,p2, points(a(nn),:), points(b(nn),:),'EdgeCross' );
        end        
        crossingCount = sum( sign(crossings) );
        if mod(crossingCount,2)
            % There where even number of crossings
           count = count+1;
        end
        if ~mod(n,1e3)
            waitbar(n/N,hWaitBar);
        end
        n = n+1;
    end % while n < N
    if ishandle(hWaitBar)
        delete(hWaitBar);
        totalA = dx*dy;
        area = count/N*totalA;
    else
        area = 0;
    end


% --- Executes on selection change in mode.
function mode_Callback(hObject, eventdata, handles)
    hPanels = [handles.scalePanel ; handles.areaPanel ; handles.lengthPanel ];
    set(hPanels,'visible','off');
    n = get(hObject,'value');
    set(hPanels(n),'visible','on');
    
    h = {handles.hScaleCalibP ,...
        handles.hAreaLines , handles.hAreaPatches, handles.hAreaPoints,...
        handles.lengthLines, handles.hLengthPointsPlot, handles.hLengthData};
    
    for nn =1: length(h)
        set(h{nn},'visible','off');
    end
    switch n
        case 1
                % scale calibration data
                set(h{1},'visible','on');
        case 2
        % Area measurment data
            set(h{ 2 },'visible','on');
            set(h{ 3 },'visible','on');
            set(h{ 4 },'visible','on');
        case 3 
            % length measures data
            set(h{ 5 },'visible','on');
            set(h{ 6},'visible','on');
            set(h{ 7 },'visible','on');
    end
    
    updateTitle(handles);

% hObject    handle to mode (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: contents = get(hObject,'String') returns mode contents as cell array
%        contents{get(hObject,'Value')} returns selected item from mode

    function updateTitle(handles)
        n = get(handles.mode,'value');
        temp = get(handles.totalArea,'userdata');
        switch n
            case 1
                    % scale calibration data
                    str = ['Scale = ' num2str(temp(1)) ' units/pixel'];
                    set(handles.totalArea,'string',str);
            case 2
            % Area measurment data
                str = ['Total Area is = ' num2str(temp(2)) ];
                set(handles.totalArea,'string',str);
            case 3 
                % length measures data
                str = ['Total Length is = ' num2str(temp(3)) ];
                set(handles.totalArea,'string',str);
        end

% --- Executes during object creation, after setting all properties.
function mode_CreateFcn(hObject, eventdata, handles)
% hObject    handle to mode (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end