(* ::Package:: *)

Needs["ErrorBarPlots`"]


(* ::Section:: *)
(*Triangle exponent*)


(* ::Text:: *)
(*Expected number of triangles is the following powerlaw*)


(* ::DisplayFormula:: *)
(*#triangles = const\[Cross]n^(T(\[Tau]))      where    T(\[Tau])=3/2 (3-\[Tau])*)


(* When importing from exponent-only-data file *)
gsraw=Import[NotebookDirectory[]<>"data/graphdata_exponent.m"];
gsraw=SortBy[gsraw,#[[1,1]]&]; (* Sort by n *)
averagesGrouped=GatherBy[gsraw,{#[[1,2]]&,#[[1,1]]&}];


(* averagesGrouped[[ tau index, n index, run index , {ntau, avgtri} ]] *)
nlabels=Map["n = "<>ToString[#]&,averagesGrouped[[1,All,1,1,1]]];
taulabels=Map["tau = "<>ToString[#]&,averagesGrouped[[All,1,1,1,2]]];
averagesErrorBars=Map[
{{#[[1,1,1]],Mean[#[[All,2]]]},
ErrorBar[StandardDeviation[#[[All,2]]]]
}&,averagesGrouped,{2}];


ErrorListPlot[averagesErrorBars,Joined->True,PlotMarkers->Automatic,PlotRange->All,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels]


ListLogLogPlot[averagesErrorBars[[All,All,1]],Joined->True,PlotMarkers->Automatic,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels]


(* ::Subsection:: *)
(*Fitting the log-log-plot*)


loglogdata=Log[averagesErrorBars[[All,All,1]]];
fits=Map[Fit[#,{1,logn},logn]&,loglogdata];
fitsExtra=Map[LinearModelFit[#,logn,logn]&,loglogdata];


fitsExtra[[1]]["ParameterConfidenceIntervalTable"]
fitsExtra[[1]]["BestFitParameters"]
fitsExtra[[1]]["ParameterErrors"]
fitsExtra[[1]]["ParameterConfidenceIntervals"]


Show[ListLogLogPlot[averagesErrorBars[[All,All,1]],Joined->True,PlotMarkers->Automatic,AxesLabel->{"n","\[LeftAngleBracket]triangles\[RightAngleBracket]"},PlotLegends->taulabels],Plot[fits,{logn,1,2000}]]


(* ::Subsection:: *)
(*Plot of T(\[Tau])*)


tauValues=averagesGrouped[[All,1,1,1,2]];
exponents=Transpose[{tauValues,fits[[All,2,1]]}];
Show[ListPlot[exponents,Joined->True,PlotMarkers->Automatic,AxesLabel->{"tau","T(\[Tau])"},PlotRange->{{2,3},{0,1.6}}],Plot[3/2(3-tau),{tau,2,3}]]


(* ::Subsection:: *)
(*T(\[Tau]) including error bars*)


tauValues=averagesGrouped[[All,1,1,1,2]];
exponentsErrorBars=Map[{{#[[1]],#[[2]]["BestFitParameters"][[2]]},ErrorBar[#[[2]]["ParameterConfidenceIntervals"][[2]]-#[[2]]["BestFitParameters"][[2]]]}&,
Transpose[{tauValues,fitsExtra}]];
Show[ErrorListPlot[exponentsErrorBars,Joined->True,PlotMarkers->Automatic,AxesLabel->{"tau","T(\[Tau])"},PlotRange->{{2,3},{0,1.6}}],Plot[3/2(3-tau),{tau,2,3}]]
