%%%%%%%%%%%%%%%%%%%%%%%%%%-----------simulation area-----------%%%%%%%%%%%%%%%%%%%%%%
% devloped in Iran(ardebil)                                 %
% Copyright (C) 2017(1396) by moussa hasanzadeh.           %
%All rights reserved.                                      %
%programed in mathlab 2016                                %
%gmail:mo30no@gmail.com
%phone: 09147082079


function vlia1
        clc
        clear;
        delete(allchild(0));
        sideN=4;
        t1=linspace(0,1,sideN+1);%vector between 0 and 1
        X=-sin(2*pi*t1+(pi/sideN));%node X coordinates with phaze shift
        Y=-cos(2*pi*t1+(pi/sideN));%node Y coordinates with phaze shift
        X=((((X-min(X))/(max(X)- min(X)))*2)-1)./5;%normalize  and schale X
        Y=((((Y-min(Y))/(max(Y)- min(Y)))*2)-1)./5;%normalize and schale Y
        Z= (1:length(Y))*0 ;

        for x=0:20:40
            for y=0:20:40
                for z=1:15
                    plot3(X+x,Y+y,Z+z,'k');hold on
                    plot3([X+x;X+x],[Y+y;Y+y],[Z+z*0;Z+z],'k');hold on
                end
            end
        end



        xVec=[0,20,40,40,20,0,0];
        yVec=[0,0,0,40,40,40,0];
        zVec1=yVec*0;
        zVec1([2,5])=zVec1([2,5])+3;
        [xVec,yVec]=dividerE(xVec,yVec,10);

        zVec=yVec*0;
        [~,zVec1]=dividerE(xVec,zVec1,10);



        hj1=fill3([xVec],[yVec],[zVec+5],'w');alpha(hj1 ,1)
        hj1=fill3([xVec],[yVec],[zVec+5+0.2],'k');alpha(hj1 ,0.2)

        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+5;zVec+5+0.2],'k');alpha(hj3 ,0.7)

        hj2=fill3([xVec],[yVec],[zVec+10],'w');alpha(hj2 ,1)
        hj2=fill3([xVec],[yVec],[zVec+10+0.2],'k');alpha(hj2 ,0.2)

        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+10;zVec+10+0.2],'k');alpha(hj3 ,0.7)

        hj3=fill3([xVec],[yVec],[zVec+15],'k');alpha(hj3 ,0.2)
        hj3=fill3([xVec;xVec],[yVec;yVec],[zVec+15;zVec+15+0.2],'k');alpha(hj3 ,0.7)

        xVecPentHous=[0,20,40,20,0];
        yVecPentHous=[0,0,0,0,0] ;
        zVecPentHous=[15,15,15,18,15];
        for mm=1:39
            % fill3(xVecPentHous,yVecPentHous+mm,zVecPentHous,'k--')
            km=fill3(xVecPentHous,yVecPentHous+mm,zVecPentHous,'g-');hold on
            alpha(km,0.2)
        end

        fill3([xVec],[yVec],[zVec+15],'r')



        zVec1=zVec1+15;

        fill3([xVec],[yVec],[zVec1],'r-','edgecolor','k')
        fill3([xVec],[yVec],[zVec1-0.2],'r-','edgecolor','k')

        hj1=fill3([xVec*2]-20,[yVec*2]-20,[zVec],'g');alpha(hj1 ,1)
        yardVect=[13:length(xVec),1:10]

        hj1=fill3([[xVec(yardVect)*2]-20;[xVec(yardVect)*2]-20],[[yVec(yardVect)*2]-20;[yVec(yardVect)*2]-20],[zVec(yardVect);zVec(yardVect)+2],'r');




        for difIndex=1:10
            zVec=zVec+0.2;
            hj1=plot3([xVec(yardVect)*2]-20,[yVec(yardVect)*2]-20,[zVec(yardVect)],'ks-','lineWidth',1);alpha(hj1 ,1);hold on

        end
        zVec(11:12)=zVec(11:12)+1;
        xVec(11)= xVec(10);
        xVec(12)= xVec(13);
       
        hj1=plot3([xVec*2]-20,[yVec*2]-20,[zVec],'k-','lineWidth',2);alpha(hj1 ,1);hold on
         zVec([10,13])=0;
        
         door=[10,11,12,13,10]
         matrixObj.x=[xVec(door)*2]-20;
         matrixObj.y=[yVec(door)*2]-20;
         matrixObj.z=[zVec(door)];
         [ newx, newy, newz]=rotateParameter(matrixObj,[0,0,1],90);
         
         hj1=fill3(newx-sqrt((xVec(10)-xVec(13)).^2)*1.3, newy-32, newz ,'w','lineWidth',2);alpha(hj1 ,1);hold on

        

        end


        function   [X,Y]=dividerE(X,Y,n)

        T1=0;
        T2=0;

        for i=1:length(X)-1
            temp1=0;
            temp1=linspace(X(i),X(i+1),n+2)  ;

            T1=[T1,temp1(1:end-1)];
        end
        T1=[T1,temp1(end)];
        X=T1(2:end);


        for i=1:length(Y)-1
            temp2=0;
            temp2=linspace(Y(i),Y(i+1),n+2)  ;

            T2=[T2,temp2(1:end-1)];
        end
        T2=[T2,temp2(end)];
        Y=T2(2:end);

        end






            function [ newx, newy, newz]=rotateParameter(matrixObj,azel,alpha)

            u = azel(:)/norm(azel);
            alph = alpha*pi/180;
            cosa = cos(alph);
            sina = sin(alph);
            vera = 1 - cosa;
            x = u(1);
            y = u(2);
            z = u(3);
            rot = [cosa+x^2*vera x*y*vera-z*sina x*z*vera+y*sina; ...
                x*y*vera+z*sina cosa+y^2*vera y*z*vera-x*sina; ...
                x*z*vera-y*sina y*z*vera+x*sina cosa+z^2*vera]';

            x = matrixObj.x;
            y = matrixObj.y;
            z = matrixObj.z;

            [m,n] = size(z);
            if numel(x) < m*n
                [x,y] = meshgrid(x,y);
            end

            [m,n] = size(x);
            newxyz = [x(:) , y(:) , z(:) ];
            newxyz = newxyz*rot;
            newx =   reshape(newxyz(:,1),m,n);
            newy =  reshape(newxyz(:,2),m,n);
            newz =  reshape(newxyz(:,3),m,n);

            end

این هم یک تصویر از نتیجه کار: