function [delka_zaznamu,Prumerna_vychylka_x_abs ,Prumerna_vychylka_y_abs,Prumerna_vychylka_cop_abs,delka_pohybu_cop ,avg_rychlost_CoP,delka_hl_osy , delka_vedl_osy,obsah,RMS ] = data_zebris (b)      
%close all ; clear all; clc ; b=NacitaniDat('s11m4b.xml');  
delka=length(b);
x=b;
cm = 0.8469 ;  mm = cm * 10;        % upravení aby byl graf v mm
z=10;    % zvětšení 
fvz=92;
close all 
clc
delka_zaznamu= length(b)/fvz
velikost_os_v_grafu=30; 
g=velikost_os_v_grafu;
kk=0;
%% vycentrování všech matic 


[rc,cc] = ndgrid(1:size(b{1,1},1),1:size(b{1,1},2));    %najde CoM 
    Mt = sum(b{1,1}(:));
    osa_y(1) = sum(b{1,1}(:) .* rc(:)) / Mt;
    osa_x (1)= sum(b{1,1}(:) .* cc(:)) / Mt;
   [radky_matice(1),sloupce_matice(1)]=size(b{1,1});
   [rc,cc] = ndgrid(1:size(b{1,2},1)  ,  1:size (b{1,2},2));    %najde CoM 
    Mt = sum(b{1,2}(:));
    osa_y(2) = sum(b{1,2}(:) .* rc(:)) / Mt;
    osa_x (2)= sum(b{1,2}(:) .* cc(:)) / Mt;
   [radky_matice(1),sloupce_matice(1)]=size(b{1,1});
for j=2:delka
           [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2)); 
             Mt = sum(b{1,j}(:));
             osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
             osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
             [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
        if abs(osa_y(j)-osa_y(j-1)) >  0.5
            b{1,j}= [zeros(1,sloupce_matice(j)) ; b{1,j}];         % přidá nulový řádek 
            [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));     
            Mt = sum(b{1,j}(:));
            osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
            [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
             if abs(osa_y(j)-osa_y(j-1)) >  0.5
                b{1,j}= [zeros(1,sloupce_matice(j)) ; b{1,j}];
                [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));     
                Mt = sum(b{1,j}(:));
                osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
                [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
                  if abs(osa_y(j)-osa_y(j-1)) >  0.5
                       b{1,j}= [zeros(1,sloupce_matice(j)) ; b{1,j}];
                       [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));     
                       Mt = sum(b{1,j}(:));
                      osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
                      [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
                     if abs(osa_y(j)-osa_y(j-1)) >  0.5
                      
                           b{1,j}= [zeros(1,sloupce_matice(j)) ; b{1,j}];
                            [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));     
                             Mt = sum(b{1,j}(:));
                             osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
                             [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
%                          
                          if abs(osa_y(j)-osa_y(j-1)) >  0.5
                                    b{1,j}(1:4,:)=[];
                                    [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));     
                                    Mt = sum(b{1,j}(:));
                                   
                                    osa_y(j) = sum(b{1,j}(:) .* rc(:)) / Mt;
                                    [radky_matice(j),sloupce_matice(j)]=size(b{1,j});
                         end; 
                       end ;  
                  end ;   
             end ;  
        end;

 [radky_matice(j),sloupce_matice(j)]=size(b{1,j}); 
if abs(osa_x(j)-osa_x(j-1)) >=  abs(0.5)
    b{1,j}= [zeros(radky_matice(j),1)   b{1,j}]  ;            
    [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));
    Mt = sum(b{1,j}(:));
    osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
     if abs(osa_x(j)-osa_x(j-1)) >=  abs(0.5)
        b{1,j}= [zeros(radky_matice(j),1)   b{1,j}]  ;            
         [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));
         Mt = sum(b{1,j}(:));
          osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
        if abs(osa_x(j)-osa_x(j-1)) >=  abs (0.5)
            b{1,j}= [zeros(radky_matice(j),1)   b{1,j}]  ;            
             [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));
             Mt = sum(b{1,j}(:));
              osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
                if abs(osa_x(j)-osa_x(j-1)) >=  abs (0.5)
                    b{1,j}= [zeros(radky_matice(j),1)   b{1,j}]  ;            
                    [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));
                     Mt = sum(b{1,j}(:));
                     osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
                    if abs(osa_x(j)-osa_x(j-1)) >=  abs (0.5)
                        b{1,j}(:,1:4)=[]  ;             
                        [rc,cc] = ndgrid(1:size(b{1,j},1)  ,  1:size (b{1,j},2));
                         Mt = sum(b{1,j}(:));
                         osa_x (j)= sum(b{1,j}(:) .* cc(:)) / Mt;
                          
                    end;  
                end; 
          end;  
     end; 
 end;
end

prum_x = mean(osa_x);             % pruměrná nulová výchylka na X
prum_y = mean(osa_y) ;            % pruměrná nulová výchylka na Y
stred = [prum_x,prum_y]           % stred v matici
nulova_vychylka_x = (osa_x - prum_x) * mm;    % nulová vychylka pro graf jednotka je mm
nulova_vychylka_y = (osa_y - prum_y) * mm;
zvetseny=imresize((x{1,520}), z);  % pro vykreslení zvětšené matice 

%% Výchylka CoP

vychylka_cop = sqrt ((nulova_vychylka_x).^2 + (nulova_vychylka_y).^2);
RMS=sqrt((1/delka)*sum((vychylka_cop.^2)));
% číselná hodnota pruměrné výchylky která je v absolutní hodnotě % v mm 
Prumerna_vychylka_x_abs    = mean (abs(nulova_vychylka_x)   )      
Prumerna_vychylka_y_abs    = mean (abs(nulova_vychylka_y) )
Prumerna_vychylka_cop_abs  = mean (vychylka_cop)


krok=delka_zaznamu/delka;
t=0.00000001:krok:delka_zaznamu;

%%  vypočet délky a rychlosti CoP 
nulova_vychylka_xx=medfilt1(nulova_vychylka_x,19);
nulova_vychylka_yy=medfilt1(nulova_vychylka_y,19);
dlk=0;
for i=1:delka-1
   cv = sqrt((nulova_vychylka_xx(1,i)- nulova_vychylka_xx(1,i+1))^2 + (nulova_vychylka_yy(1,i)- nulova_vychylka_yy(1,i+1))^2);  % vzdalenost bodu i a i+1
   rychlost(i)=cv;   
   dlk=dlk+cv;    % suma cv  
end
delka_pohybu_cop=dlk   
   avg_rychlost_CoP=mean(rychlost)*fvz  % průměpruměxrná rychlost COP mm

   %% Výpočet komf elipsy 95%
  % v případě že chci vykreslit konfidenční elipsu do rozožení tlaku musím použít
%  osa_x, osa_y a potom musím na řádku 177 a 178 odpoznámkovat konstantu mm
y1=nulova_vychylka_x';             
y2=nulova_vychylka_y';
% y1=osa_x';
% y2=osa_y';
data = [y1 y2];
kovariance = cov(data);
[vlastni_vektor, vlastni_cislo ] = eig(kovariance);

[nejvetsi_vlast_vekt_radek, s] = find(vlastni_cislo == max(max(vlastni_cislo)));  
nejvetsi_vlast_vek = vlastni_vektor(:, nejvetsi_vlast_vekt_radek);

nejvetsi_vlast_cislo = max(max(vlastni_cislo));
    if(nejvetsi_vlast_vekt_radek == 1)
      nejmensi_vlast_cislo = max(vlastni_cislo(:,2));
      nejmensi_vlast_vekt = vlastni_vektor(:,2);
    else
      nejmensi_vlast_cislo = max(vlastni_cislo(:,1));
      nejmensi_vlast_vekt = vlastni_vektor(1,:);
    end

uhel = atan2(nejvetsi_vlast_vek(2), nejvetsi_vlast_vek(1));  %Vypočte úhel mezi osou X a nejvetsim vlastním vektorem
    if(uhel < 0)
      uhel = uhel + 2*pi;    %posune sinusovku z (-pi až pi) do (0 až 2 pi)
    end
prumer = mean(data);
chi_kvadrat = 2.4477;             
sinusovka = linspace(0,2*pi);
X0=prumer(1);              % střed elipsy
Y0=prumer(2);              % střed elipsy
a=chi_kvadrat*sqrt(nejvetsi_vlast_cislo); % standardní odchylka (dosadíme li za velikost poloos její dvojnásobek dostaneme 95%comf el/ dosadíme li 3 násobek máme 99% comf el
c=chi_kvadrat*sqrt(nejmensi_vlast_cislo);
elipsa_x  = a*cos( sinusovka );
elipsa_y  = c*sin( sinusovka );
R = [ cos(uhel) sin(uhel); -sin(uhel) cos(uhel) ];  % Definuje správné natočení elipsy
elipsa = [elipsa_x;elipsa_y]' * R;
%%výoočet délky poloos 
delka_hl_osy = 2*a%*mm   % v případě že dám osa_x/y 
delka_vedl_osy = 2*c%*mm
obsah = pi * (delka_hl_osy/2) * (delka_vedl_osy/2)
  if delka_hl_osy > delka_vedl_osy   
               delka_hl_osy = delka_hl_osy;
               delka_vedl_osy = delka_vedl_osy;
           else
               delka_hl_osy = delka_vedl_osy;
               delka_vedl_osy = delka_hl_osy;
           end
              




%% GRAFY 

% 
% imshow(zvetseny);
%  
%  figure
%  plot(t,vychylka_cop )   
% axis([0 delka_zaznamu 0 g])                  
% title('Stabilogram') 
%    xlabel('t [s]')                            
%    ylabel('vychylka osy x [mm]') 
%    figure                          
% plot(t,nulova_vychylka_x )
% axis([0 delka_zaznamu -g g])
% title('Pohyb CoP v ose x v závislosti na čase') 
%    xlabel('t [s]') 
%    ylabel('vychylka osy x [mm]')  
% figure
% plot(t,nulova_vychylka_y )  
% axis([0 delka_zaznamu -g g])
% title('Pohyb CoP v ose y v závislosti na čase') 
%    xlabel('t [s]') 
%    ylabel('vychylka osy y [mm]')  
%    
%    figure
% plot(t,vychylka_cop )  
% axis([0 delka_zaznamu 0 g])
% title('Stabilogram') 
%    xlabel('t [s]') 
%    ylabel('vychylka CoP [mm]')  
%      
% figure 
% plot(nulova_vychylka_x(1,:) , nulova_vychylka_y(1,:))
% axis([-g g -g g])
% title('Statokinesigram') 
%    xlabel('osa X [mm]')
%    ylabel('osa Y [mm]' )
%    
%    figure
%  plot(elipsa(:,1) + X0,elipsa(:,2) + Y0,'-')
% hold on;
%  plot(medfilt1(data(:,1)), medfilt1(data(:,2)) );
% axis equal
% 
% hold on;
% hXLabel = xlabel('x');
% hYLabel = ylabel('y');
% 
% hXLabel = xlabel('x');
% hYLabel = ylabel('y');
%    
% 
% 
% figure    
% image(zvetseny*5)
%   hold on
% plot(data(:,1)*z, data(:,2)*z,'rx'); 
% hold on 
% plot((elipsa(:,1) + X0)*z,(elipsa(:,2) + Y0)*z,'y-')
end