SetDirectory[NotebookDirectory[]];
data=Import["D01_2-1x16-25bbb.xls"]; 
(* importovani dat od Cez Distribucni sluzby s.r.o ze souboru *)

pocetlidi=150;
 (* celkovy pocet koncovych odberatelu v namerenych datech *)
proud=16;
(* velikost proudoveho jistice *)
napeti=230;
(* napeti *)

NactiSloupec[i_]:=Table[data[[1,k,i]],{k,2,2976}];
 (*funkce, ktera nacita sloupec*)
(*nacteni hodnot*)
GrafClovek[k_,proud_,napeti_]:=Module[{clovek3,b}, (* GrafClovek je funkce, vybere k-teho odberatele a prevede kWh na pomerne hodnoty  *)
clovek3=NactiSloupec[k];(* vstupem GrafClovek je k(ktery clovek) proud a napeti *)

b=Table[(1/4)(clovek3[[4*i+1]]+clovek3[[4*i+2]]+clovek3[[4*i+3]]+clovek3[[4*i+4]])-(1/4)(clovek3[[4*(i-1)+1]]+clovek3[[4*(i-1)+2]]+clovek3[[4*(i-1)+3]]+clovek3[[4*(i-1)+4]]),{i,1,Round[(2975/4)]-2,1}];(*b - prepocet z 15 minutoveho mereni na hodinove mereni*)

((proud*napeti)/1000)*b (*prepocet na pomerne hodnoty*)

]

(* ListPlot[GrafClovek[3,16,230],Joined->True] (* test funkce GrafClovek *)

diagramdata=Import["Normalizovane_TDD_2013_CZ.xls"];(*import referenčcnich TDD z ote-cr.cz*)

diagramDataGraf=Table[diagramdata[[1,i,7]],{i,7,7+741}];(*vybrani sloupce s TDD4 z ote-cr.cz *)

grafTDD4=ListPlot[diagramDataGraf,Joined->True,PlotStyle->Red,PlotLabel->"TDD4",AxesLabel->{"Hodina [h]","Pomerna hodnota [-]"}](*zobrazi referencni graf TDD4*)

nahData=Table[diagramDataGraf+0.1Random[Real,{-1,1}],{i,1,150}];
 (* testovaci data, ktera porovnava referencni data s referencnimi daty. Pomoci fce Random se o nahodnou hodnotu z intervalu [-1,1] zasumime referencni data. Pro testovani se musi zamenit GrafClovek za GrafClovek2 (data z  www.ote-cr.cz, TDD4 ) *)
GrafClovek2[k_,proud,napeti]:=nahData[[k]]
fcePrumerProTDD[k_,n_]:=Module[{poradiList,plusy,vyst,listy,listySoucet,soucetN},
(* nazev funkce, ktera zpracovava data a vytvori sloupec pro porovnani s referencnimi daty (data z  www.ote-cr.cz, TDD4 ) *)
poradiList=Table[Flatten[RandomPermutation[pocetlidi][[1]]],{i,1,k}];
(*vytvoreni nahodne permutace - zameneni sloupcu a nasledny vyber od prvniho*)
listy=Table[Table[GrafClovek[poradiList[[b,i]],proud,napeti],{i,1,n}],{b,1,k}];
(* vytvoreni k kombinaci z poctu konecnych odberatelu; testovaci program - zamenit z GrafClovek na GrafClovek2  *)
listySoucet=((Plus@@Table[listy[[i]],{i,1,k}])/k);
(* prumer k kombinaci poctu konecnych odberatelu(n) vytvorenych v predchozim kroku *)
soucetN=((Plus@@Table[listySoucet[[a]],{a,1,n}])/n)
(* prumer konecnych odberatelu - vytvoreni jednoho sloupce *)

]

fceRozdilTDDaData[k_,n_]:=Module[{},
(Sqrt[((Plus@@(((fcePrumerProTDD[k,n]-diagramDataGraf))^2))/742)])
]
(* vypocet stredni kvadraticke odchylky *)
datavysledek =Table[Table[fceRozdilTDDaData[k,n],{n,1,140}],{k,5,85,20}];

(* vysledky pro k=5, k=25, k=45, k=65, k=85; konstanta k je pocet permutaci sloupcu v namerenych datech od CEZ Distribucni sluzby s.r.o. (konecnych odberatelu) *)
(* pro zkousku je nutne zrusit vnejsi Table; pro nizky cas vypoctu je zapotřebí zvolit k = 1 *)

Joined -> True, AxesLabel{"n [Pocet lidi]", "stredni kvadraticka odchylka [pomerne jednotky]"}]

porovnaniTDDkje5=ListPlot[datavysledek[[1]],AxesLabel->{"n [Pocet koncovych odbernych mist]", "stredni kvadraticka odchylka [pomerne jednotky]"},PlotLabel->"k = 5",Joined->True]
porovnaniTDDkje25=ListPlot[datavysledek[[2]],AxesLabel->{"n [Pocet koncovych odbernych mist]", "stredni kvadraticka odchylka [pomerne jednotky]"},PlotLabel->"k = 25",Joined->True]
porovnaniTDDkje45=ListPlot[datavysledek[[3]],AxesLabel->{"n [Pocet koncovych odbernych mist]", "stredni kvadraticka odchylka [pomerne jednotky]"},PlotLabel->"k = 45",Joined->True]
porovnaniTDDkje65=ListPlot[datavysledek[[4]],AxesLabel->{"n [Pocet koncovych odbernych mist]", "stredni kvadraticka odchylka [pomerne jednotky]"},PlotLabel->"k = 65",Joined->True]
porovnaniTDDkje85=ListPlot[datavysledek[[5]],AxesLabel->{"n [Pocet koncovych odbernych mist]", "stredni kvadraticka odchylka [pomerne jednotky]"},PlotLabel->"k = 85",Joined->True]

fceVytvorTDD[k_,n_]:=Module[{},
(* funkce, ktera vytvori TDD pro dane k a n, kde k pocet kombinaci n koncovych odberatelu. k je pocet cyklu a n je pocet koncovych odberatelu*)
ListPlot[fcePrumerProTDD[k,n],Joined->True,PlotLabel->"vytvorene TDD pro n=50, k = 85",
AxesLabel->{"Hodina [h]", "Pomerna hodnota[-]"}]
]

fceOdchylkaTDD[k_,n_]:=Module[{nty,nMinusPrvni,rozdil,kvadrat,soucetKvadratuOdchylek,odmocnina,podil},
(* funkce, ktera vytvori zavislost procentualni odchylky pro n lidi a k permutaci *)
nty=fcePrumerProTDD[k,n];
(* definovani n koncoveho odberatele *)
nMinusPrvni=fcePrumerProTDD[k,n-1];
(* definovani n-1 odberatele *)
rozdil=nty-nMinusPrvni;
(* rozdil dat n koncoveho odberatele a n-1 koncoveho odberatele *)
kvadrat=(rozdil^2);
(* umocneni funkce rozdil, zpresneni odchylky *)
soucetKvadratuOdchylek=Plus@@kvadrat;
(* soucet kvadratyckych odchylek *)
odmocnina=(Sqrt[soucetKvadratuOdchylek]);
(* odmocneni  *)

podil=(odmocnina/n)*100
(* prumerna smerodatna odchylka vztazena na jednoho odberatele *)
]

dataLidi=Table[Table[{i,fceOdchylkaTDD[k,i]},{i,2,140,1}],{k,5,85,20}]
(* vypocet pro k=5,25,45,65,85 a n=140; konstanta k je pocet permutaci sloupcu v namerenych datech od CEZ Distribucni sluzby s.r.o. (konecnych odberatelu)  *)

vytvoreneTDD=fceVytvorTDD[85,50]
(* vytvoreni grafu TDD pro hodnoty k=85 a n=50 *)

kje5=ListPlot[dataLidi[[1]],AxesLabel->{"n [pocet koncovych odberatelu]", "prumerna smerodatna odchylka na jednoho odberatele [%]"},PlotLabel->"k = 5",Joined->True]
kje25=ListPlot[dataLidi[[2]],AxesLabel->{"n [pocet koncovych odberatelu]", "prumerna smerodatna odchylka na jednoho odberatele [%]"},PlotLabel->"k = 25",Joined->True]
kje45=ListPlot[dataLidi[[3]],AxesLabel->{"n [pocet koncovych odberatelu]", "prumerna smerodatna odchylka na jednoho odberatele [%]"},PlotLabel->"k = 45",Joined->True]
kje65=ListPlot[dataLidi[[4]],AxesLabel->{"n [pocet koncovych odberatelu]", "prumerna smerodatna odchylka na jednoho odberatele [%]"},PlotLabel->"k = 65",Joined->True]
kje85=ListPlot[dataLidi[[5]],AxesLabel->{"n [pocet koncovych odberatelu]", "prumerna smerodatna odchylka na jednoho odberatele [%]"},PlotLabel->"k = 85",Joined->True]

Export["kje5.tiff",kje5]
Export["kje25.tiff",kje5]
Export["kje45.tiff",kje5]
Export["kje65.tiff",kje5]
Export["kje75.tiff",kje5]
Export["vytvoreneTDD.tiff",vytvoreneTDD]
Export["porovnaniTDDkje5.tiff",porovnaniTDDkje5]
Export["porovnaniTDDkje25.tiff",porovnaniTDDkje25]
Export["porovnaniTDDkje45.tiff",porovnaniTDDkje45]
Export["porovnaniTDDkje65.tiff",porovnaniTDDkje65]
Export["porovnaniTDDkje85.tiff",porovnaniTDDkje85]

(* export vytvorenych grafu *)