(* ::Package:: *)

Needs["ErrorBarPlots`"]


(* ::Section:: *)
(*Data import*)


gsraw=Import[NotebookDirectory[]<>"data/graphdata_dsp.m"];
(* gsraw=SortBy[gsraw,{#[[1,1]]&,#[[1,2]]&}]; (* Sort by n and then by tau. The {} forces a *stable* sort because otherwise Mathematica sorts also on triangle count and other things. *) *)


gdata=GatherBy[gsraw,{#[[1,2]]&,#[[1,1]]&}];
(* Data format: *)
(* gdata[[ tau index, n index, run index , datatype index ]] *)
(* datatype index:
1: {n,tau}
2: avgtris
3: dsp
*)
nlabels=Map["n = "<>ToString[#]&,gdata[[1,All,1,1,1]]];
taulabels=Map["\[Tau] = "<>ToString[#]&,gdata[[All,1,1,1,2]]];


(* ::Section:: *)
(*Plot triangle counts vs dsp*)


(* ::Subsection:: *)
(*Compute dsp from ds*)


getProperty[ds1_]:=Module[{ds,n=Length[ds1],tmp=ConstantArray[0,{Length[ds1],Length[ds1]}]},
ds=N[ds1/Sqrt[N[Total[ds1]]]]; (* scale degrees by 1/Sqrt[total] *)
(* The next table contains unneeded entries, but its faster to have a square table for the sum *)
tmp=Table[1.-Exp[-ds[[i]]ds[[j]]],{i,1,n},{j,1,n}];
Sum[tmp[[i,j]]*tmp[[j,k]]*tmp[[i,k]],{i,3,n},{j,2,i-1},{k,1,j-1}] (* somehow i>j>k is about 60x faster than doing i<j<k !!! *)
(* This sparser table is slower
tmp=Table[1.-Exp[-ds[[i]]ds[[j]]],{i,1,n-1},{j,i+1,n}];
(* tmp[[a,b]] is now with  ds[[a]]*ds[[a+b]] *)
Sum[tmp[[i,j-i]]*tmp[[j,k-j]]*tmp[[i,k-i]],{i,1,n-2},{j,i+1,n-1},{k,j+1,n}]
*)
];


(* gdata[[ tau index, n index, run index , {ntau, #tris, ds} ]] *)
avgAndProp=ParallelMap[{getProperty[#[[3]]],Mean[#[[2,1;;-1]]]}&,gdata[[2,2,1;;100]]];


(* ::Subsection:: *)
(*Plot #trianges vs some degree-sequence-property*)


allPlots=Map[Show[
ListPlot[#[[All,{2,3}]],
AxesOrigin->{0,0},
Frame->True,FrameLabel->{"average number of triangles","degree-sequence-property"},
AspectRatio->Automatic,
ImageSize->Automatic,
PlotLabel->"n = "<>ToString[#[[1,1,1]]]<>" , \[Tau] = "<>ToString[#[[1,1,2]]],
PlotStyle->Black
],Plot[x,{x,1,10000}]]
&,gdata,{2}]


combiPlot=Show[
ListLogLogPlot[gdata[[All,1,All,{3,2}]],
(*AxesOrigin->{0,0},*)
(*PlotRange\[Rule]{{0,3000},{0,3000}},*)
PlotRange->Automatic,
Frame->True,FrameLabel->{"degree-sequence-property","average number of triangles"},
AspectRatio->Automatic,
ImageSize->Automatic,
PlotLabel->"n = "<>ToString[gdata[[1,1]][[1,1,1]]],
PlotLegends->taulabels
],Plot[x,{x,1,10000},PlotStyle->Black]]


plotgrid=GraphicsRow[{GraphicsColumn[allPlots[[{1,2},1]]],allPlots[[3,1]]},Spacings->0,ImageSize->600]


Export[NotebookDirectory[]<>"plots/dsp_correlations.pdf",combiPlot]