function Querp
%QUERP für Matlab
%
%Ermittlung der Querschnittswerte von Flächen mit abschnittsweise
%geradliniger Begrenzung mit Hilfe elektronischer Rechenanlagen
%
%Autor: Gerhard Knappstein
%Matlab-Bearbeitung: Philipp Fleiger
%
%1987 - 2007
%


men_contr=0;    %Menue-Kontrollvariable
beenden=0;

Punkte=[];
einheit=['m'];


%--------------------------------------------------------------------------
%Textmenue
%--------------------------------------------------------------------------

while beenden == 0

    clc
    
    disp('----------------------------------------------------------------- ');
    disp('**************************');
    disp('*  Q U E R P für Matlab  *');
    disp('**************************');
    disp(' ');
    disp('Universität Siegen    Lehrstuhl für numerische Mechanik ');
    disp('G. Knappstein, P. Fleiger  2007 ');
    disp('----------------------------------------------------------------- ');
    disp(' ');
    disp(' ');

    switch men_contr
        
        case 0
            disp('( 1 ) Eckpunkte eingeben');
            disp('( 2 ) Graphische Eingabe');
            disp('( 3 ) Beenden');
            
        case 1
            disp('( 1 ) Weitern Punkt hinzufügen');
            disp('( 2 ) Einen bestehenden Punkt ändern');
            disp('( 3 ) Ausgabe der Eckpunkte');
            disp('( 4 ) Berechnung starten (Eingabe der Punkte abschließen)');
            disp('( 5 ) Beenden');
            
        case 2
            disp('( 1 ) Ausgabe der Berechnungsdaten');
            disp('( 2 ) Ausgabedatei erzeugen');
            disp('( 3 ) Punkte ändern/ hinzufügen (Zurück zur Punkteingabe)');
            disp('( 4 ) Graphische Ausgabe');
            disp('( 5 ) Einheit wählen');
            disp('( 6 ) Beenden');           
        
    end  %switch
    
    disp(' ');
    disp(' ');
    Auswahl=input('Bitte wählen Sie eine Option aus:   ');
    
    switch men_contr
        
        case 0
            
            switch Auswahl
                
                case 1
                    clc
                    Punkte=punkt_eingeben(Punkte);
                    men_contr=1;
                    
                case 2
                    clc
                    Punkte=graph_ein;
                    men_contr=1;
                    
                case 3
                    beenden=1;
                    
                otherwise
                    disp(' ');
                    disp('Ungültige Eingabe! Bitte versuchen Sie es erneut.')
                    pause
            end
            
        case 1
            
            switch Auswahl
                
                case 1
                    clc
                    Punkte=punkt_eingeben(Punkte);
                    
                case 2
                    clc
                    Punkte=punkt_aendern(Punkte);
                    
                case 3
                    clc
                    ausgeben(Punkte);                    
                  
                case 4
                    clc
                    Flaeche=berechnen(Punkte);
                    men_contr=2;                
                    
                case 5
                    clc
                    beenden=1;
                    
                otherwise
                    disp(' ');
                    disp('Ungültige Eingabe! Bitte Taste drücken und erneut versuchen.')
                    pause
            end
            
        case 2
            
            switch Auswahl
                
                case 1
                    clc
                    ausgabe_ges(Flaeche, Punkte, einheit);
                    
                case 2
                    clc
                    datei(Flaeche, Punkte);
                    
                case 3
                    clc
                    men_contr=1;
                    
                case 4
                    graph_aus(Flaeche, Punkte);
                    
                case 5
                    clc
                    einheit=ein_waehlen;
                    
                case 6
                    beenden=1;
                    
                otherwise
                    disp(' ');
                    disp('Ungültige Eingabe! Bitte Taste drücken und erneut versuchen.')
                    pause
            end
                    
        
    end  %switch men_contr
           



end %Menue


%--------------------------------------------------------------------------
%Unterfunktionen
%--------------------------------------------------------------------------



    function Punkte=punkt_eingeben(Punkte)
        clc
        
        disp(' ');
        disp('Bitte geben Sie einen Eckpunkt an!')
        disp(' ');
        Py=input('y_quer-Koordinate:  ');
        Pz=input('z_quer-Koordinate:  ');
        Punkte=[Punkte; Py, Pz];
        pause
    end

    function Punkte=graph_ein
        clc

        %Zeichenfläche erzeugen

        disp('Bitte ungefähre Größe der Fläche angeben.');
        disp(' ');
        ymax=input('y-max: ');
        zmax=input('z-max: ');

        figure
        axis([0, ymax, 0, zmax]);
        hold on
        %axis equal


        %Zeichnen

        [y, z, Taste]=ginput(1);
        Punkte=[y, z];

        while Taste ~= 3
            [y, z, Taste]=ginput(1);
            y_draw=[Punkte(end,1), y];
            z_draw=[Punkte(end,2), z];
            Punkte=[Punkte; y, z];
            line(y_draw, z_draw);
        end

        y_draw=[Punkte(end,1), Punkte(1,1)];
        z_draw=[Punkte(end,2), Punkte(1,2)];
        line(y_draw, z_draw);
        
        pause

    end

    function Punkte=punkt_aendern(Punkte)
        ausgeben(Punkte);
        disp(' ');
        disp(' ');
        disp('Welcher Punkt soll geändert werden? ');
        i=input('Nummer der Punktes:   ');
        disp(' ');
        disp('Bitte geben Sie die neuen Koordinaten an!')
        disp(' ');
        Py=input('y_quer-Koordinate:  ');
        Pz=input('z_quer-Koordinate:  ');
        Punkte(i,:)=[Py, Pz];
        
    end

    function ausgeben(Punkte)
        clc
        disp(' ');
        disp('Bisher eingegebene Punkte: ');
        disp(' ');
        for i=1:length(Punkte)
            ausg=[num2str(i), '   ', num2str(Punkte(i,1)), '   ', num2str(Punkte(i,2))];
            disp(ausg);
        end
        pause
    end

    function Flaeche=berechnen(Punkte)
        A=polyarea(Punkte(:,1), Punkte(:,2));   %Berechnung der Fläche
        Punkte=[Punkte; Punkte(1,:)];            %Anfangspunkt=Endpunkt
        
        %Statische Momente
        
        Sq_y=0;
        Sq_z=0;
        
        for i=1:length(Punkte)-1
            Sq_y=Sq_y+1/6*((Punkte(i,1)*Punkte(i+1,2)-Punkte(i+1,1)*Punkte(i,2))*(Punkte(i,2)+Punkte(i+1,2)));
            
            Sq_z=Sq_z+1/6*((Punkte(i,1)*Punkte(i+1,2)-Punkte(i+1,1)*Punkte(i,2))*(Punkte(i,1)+Punkte(i+1,1)));
        end
        
        %Koordinaten des Schwerpunktes
        
        zq_s=Sq_y/A;
        yq_s=Sq_z/A;
        
        %Flächenmomente 2. Ordnung in Bezug auf die Bezugsachsen y_q und z_q
        
        I_yq=0;
        I_zq=0;
        I_yzq=0;
        
        for i=1:length(Punkte)-1
            %axiale Flächenträgheitsmomente
            
            I_yq=I_yq+(1/12*((Punkte(i,1)*Punkte(i+1,2)-Punkte(i+1,1)*Punkte(i,2))*((Punkte(i,2)+Punkte(i+1,2))^2-Punkte(i,2)*Punkte(i+1,2))));
            
            I_zq=I_zq+(1/12*((Punkte(i,1)*Punkte(i+1,2)-Punkte(i+1,1)*Punkte(i,2))*((Punkte(i,1)+Punkte(i+1,1))^2-Punkte(i,1)*Punkte(i+1,1))));
            
            %Zentrifugalmoment
            
            I_yzq=I_yzq-1/12*((Punkte(i,1)*Punkte(i+1,2)-Punkte(i+1,1)*Punkte(i,2))*((Punkte(i,1)+Punkte(i+1,1))*(Punkte(i,2)+Punkte(i+1,2))-1/2*(Punkte(i,1)*Punkte(i+1,2)+Punkte(i+1,1)*Punkte(i,2))));
            
        end
        
        I_y=I_yq-zq_s^2*A;
        I_z=I_zq-yq_s^2*A;
        
        I_yz=I_yzq+yq_s*zq_s*A;
        
        Flaeche=[A, Sq_y, Sq_z, yq_s, zq_s, I_y, I_z, I_yz];
        
        
        disp('Berechnung abgeschlossen. Bitte Taste drücken, um ins Ausgabemenue zu gelangen.');
        pause
        
    end

    function ausgabe_ges(Flaeche, Punkte, einheit)
        format short g
        disp(' ');
        disp('Koordinaten der Eckpunkte (y_quer,z_quer): ');
        disp(' ');
        disp(Punkte);
        disp(' ');
        disp(' ');
        disp('Berechnete Querschnittswerte: ');
        disp(' ');
        aus=['Flaecheninhalt:                   ', num2str(Flaeche(1)), ' ', einheit, '^2'];
        disp(aus);
        disp(' ');
        disp('statische Momente')
        disp(' ');
        aus=['- bezogen auf die y_quer-Achse (Sq_y):  ', num2str(Flaeche(2)), ' ', einheit, '^3'];
        disp(aus);
        aus=['- bezogen auf die z_quer-Achse (Sq_z):  ', num2str(Flaeche(3)), ' ', einheit, '^3'];
        disp(aus);
        disp(' ');
        disp('Koordinaten des Schwerpunktes: ');
        disp(' ');
        aus=['y-Koord. yq_s                   :  ', num2str(Flaeche(4)), ' ', einheit];
        disp(aus);
        aus=['z-Koord. zq_s                   :  ', num2str(Flaeche(5)), ' ', einheit];
        disp(aus);
        disp(' ');
        disp('Flächenmomente 2. Ordnung in Bezug auf die Achsen y und z ');
        disp(' ');
        disp('axiale Flächenträgheitsmomente ');
        aus=['I_y:                              ', num2str(Flaeche(6)), ' ', einheit, '^4'];
        disp(aus);
        aus=['I_z:                              ', num2str(Flaeche(7)), ' ', einheit, '^4'];
        disp(aus);
        disp(' ');
        disp('Zentrifugalmoment ');
        aus=['I_yz:                             ', num2str(Flaeche(8)), ' ', einheit, '^4'];
        disp(aus);
                
        pause
    end

    function datei(Flaeche, Punkte)
        fid = fopen('ausgabe_querp.txt','wt');
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'Koordinaten der Eckpunkte (y_quer,z_quer): \n');
        fprintf(fid,'%s \n',' ');
        for i=1:length(Punkte)
            P=[num2str(i), '  (', num2str(Punkte(i,1)), ' , ',num2str(Punkte(i,2)), ') \n' ];
            fprintf(fid,P);
        end        
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','Berechnete Querschnittswerte: ');
        fprintf(fid,'%s \n',' ');
        aus=['Flaecheninhalt:                   ', num2str(Flaeche(1)), ' ', einheit, '^2'];
        fprintf(fid,'%s \n',aus);
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','statische Momente')
        fprintf(fid,'%s \n',' ');
        aus=['- bezogen auf die y_quer-Achse (Sq_y):  ', num2str(Flaeche(2)), ' ', einheit, '^3'];
        fprintf(fid,'%s \n',aus);
        aus=['- bezogen auf die z_quer-Achse (Sq_z):  ', num2str(Flaeche(3)), ' ', einheit, '^3'];
        fprintf(fid,'%s \n',aus);
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','Koordinaten des Schnittpunktes: ');
        fprintf(fid,'%s \n',' ');
        aus=['y-Koord. yq_s                   :  ', num2str(Flaeche(4)), ' ', einheit];
        fprintf(fid,'%s \n',aus);
        aus=['z-Koord. zq_s                   :  ', num2str(Flaeche(5)), ' ', einheit];
        fprintf(fid,'%s \n',aus);
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','Flächenmomente 2. Ordnung in Bezug auf die Achsen y und z ');
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','axiale Flächenträgheitsmomente ');
        fprintf(fid,'%s \n',' ');
        aus=['I_y:                              ', num2str(Flaeche(6)), ' ', einheit, '^4'];
        fprintf(fid,'%s \n',aus);
        aus=['I_z:                              ', num2str(Flaeche(7)), ' ', einheit, '^4'];
        fprintf(fid,'%s \n',aus);
        fprintf(fid,'%s \n',' ');
        fprintf(fid,'%s \n','Zentrifugalmoment ');
        fprintf(fid,'%s \n',' ');
        aus=['I_yz:                             ', num2str(Flaeche(8)), ' ', einheit, '^4'];
        fprintf(fid,'%s \n',aus);
        fprintf(fid,'%s \n',' ');
        status=fclose(fid);
    end

    function graph_aus(Flaeche, Punkte)
%         plot(Flaeche(4), Flaeche(5), '+');
% 
%         plot(Punkte(:,1), Punkte(:,2), 'o');
    end
    
    function einheit=ein_waehlen
        clc
        disp(' ');
        disp('Bitte wählen Sie die verwendete Einheit: ');
        disp(' ');        
        disp('(1) m ');
        disp('(2) cm ');
        disp('(3) mm ');
        disp(' ');
        e=input('Bitte Zahl eingeben:  ');
        if e==1
            einheit=['m'];
        elseif e==2
            einheit=['cm'];
        elseif e==3
            einheit=['mm'];
        else
            disp(' ');
            disp('Keine gültige Eingabe!');
        end
        pause
    end


end