import os
import numpy as np
import matplotlib.pyplot as plt

def path(str):
    return os.path.join(os.path.dirname(__file__), str)

params = {'text.usetex' : True,
        'font.size' : 9,
        'axes.labelsize': 9,
        'font.family': 'serif',
        'font.serif': ['Computer Modern'],
        'lines.linewidth': 1.0,
        'lines.markersize': 4.0}
plt.rcParams.update(params)
plt.rc('text.latex')
plt.rcParams['text.latex.preamble'] = r'\usepackage{amsmath} \usepackage{amssymb} \newcommand{\Hsquare}{\text{\fboxsep=.0pt\fbox{\rule{0pt}{1ex}\rule{1ex}{0pt}}}}'


fig, (ax1, ax2) = plt.subplots(2, 1)
fig.set_size_inches(2.78, 3.4)
fig.subplots_adjust(hspace=0.1)

N = np.array([32, 64, 128, 256, 512])

# based on max(abs(I-mI))
g4    = np.array([0.5011872336272722, 0.251188643150958, 0.1778279410038923, 0.08912509381337454, 0.050118723360, 0.026607250597988])
Ig4m  = np.array([0.0005506472953252, 0.000586596778340, 0.0005976226176908, 0.00061110235090501, 0.000617123828, 0.000620769975126])
sIg4  = np.array([0.0009739116972556, 0.000297330984130, 0.0001326446135625, 6.04437438806278e-5, 2.927084e-05,   1.442193014900e-5])
maIg4 = np.array([0.0018660537751012, 0.000512053769650, 0.0002244661068709, 8.98162868739334e-5, 4.09567574e-05, 2.023768462452e-5])
Ig4   = np.array([0.0005516748327894, 0.000586897147431, 0.0005977567637132, 0.00061114310903046, 0.000617145674, 0.000620780874059])

# based on std(I)
# g4    = np.array([0.56234132519034910, 0.17782794100389230, 0.15848931924611134, 0.08912509381337454, 0.0446683592])
# Ig4m  = np.array([0.00054242537894146, 0.00059743539260809, 0.00060052502203455, 0.00061110235090501, 0.0006179633])
# sIg4  = np.array([0.00089162421524786, 0.00034951016974832, 0.00014254869693327, 6.04437438806278e-5, 3.067031e-05])
# maIg4 = np.array([0.00169756936340323, 0.00073368285673570, 0.00026796271462830, 8.98162868739334e-5, 5.851010e-05])
# Ig4   = np.array([0.00054323144629719, 0.00059787372071365, 0.00060068594034965, 0.00061114310903046, 0.0006179785])

d4    = np.array([0.19634954084936207, 0.09817477042468103, 0.04908738521234052, 0.0245436926061702, 0.01227184630, 0.00619441075076])
Id4m  = np.array([0.00059190447348925, 0.00060879920577597, 0.00061511567329138, 0.0006189368525996, 0.00062097266, 0.00062179562541])
sId4  = np.array([0.03075387217660721, 0.00339722175316666, 0.04343834594736065, 0.0425272297585360, 0.04029526010, 0.03695982861245])
maId4 = np.array([0.05781294788905573, 0.00825310584969745, 0.08460286654538818, 0.0828775465300915, 0.07839494801, 0.07195737414149])
Id4   = np.array([0.00062612421597038, 0.00060947593497851, 0.00066715622250883, 0.0005752781734968, 0.00057927775, 0.00058335328657])

IML4 = 0.0006256666923855121
IML8 = 0.0003140611503793372
IML16 = 0.0001563805017858243

g8    = np.array([1.1885022274370187, 0.6309573444801931, 0.29853826189179, 0.21134890398366465, 0.10000000000000000, 0.05623413251903491])
Ig8m  = np.array([0.0002694879083344, 0.0002889719098467, 0.00030127732579, 0.00030459935747239, 0.00030887728631451, 0.00031057351613738])
sIg8  = np.array([0.0007181495500206, 0.0002184909569620, 7.25075082094e-5, 3.22062634663118e-5, 1.49441591506481e-5, 7.25718243694343e-6])
maIg8 = np.array([0.0015987034042133, 0.0003919048255354, 0.00012330539366, 5.18567885206883e-5, 2.22271123088298e-5, 1.01067463943456e-5])
Ig8   = np.array([0.0002697268162987, 0.0002890826729251, 0.00030131447491, 0.00030461515301685, 0.00030888233241929, 0.00031057627108841])

d8    = np.array([0.387257644921, 0.1945360081622663, 0.09772372209558, 0.04908738535811, 0.02466039337234, 0.012359474334445])
Id8m  = np.array([0.000296732995, 0.0003039894367387, 0.00030867527620, 0.00030970660451, 0.00031067341899, 0.000311166225221])
sId8  = np.array([0.046939623707, 0.0350524919117734, 0.00364433891542, 0.04834744570593, 0.04729770984118, 0.045676092845318])
maId8 = np.array([0.084648495323, 0.0642091107809958, 0.00866056010099, 0.09277458220733, 0.09055307120241, 0.087104513150755])
Id8   = np.array([0.000297340382, 0.0003235083281586, 0.00030868064494, 0.00033843950535, 0.00028691210822, 0.000288173437248])



g16    = np.array([1.9952623149688793, 1.258925411794167, 0.707945784384137, 0.3349654391578276, 0.23713737056616552, 0.11220184543019636])
Ig16m  = np.array([0.0001338085293620, 0.000141884739957, 0.000148123586042, 0.0001524610523049, 0.00015361334811614, 0.00015508932815129])
sIg16  = np.array([0.0028564531715254, 0.000732676672992, 0.000126077959714, 3.3568814816081e-5, 1.36123142008695e-5, 5.05186864301191e-6])
maIg16 = np.array([0.0060330657527504, 0.001316082456436, 0.000287869470502, 5.1259233083596e-5, 2.01784890088236e-5, 7.18612870726867e-6])
Ig16   = np.array([0.0001338009209540, 0.000141966371376, 0.000148166226300, 0.0001524688673415, 0.00015361644780139, 0.00015509024995997])

d16    = np.array([0.779830110523258, 0.390840895792402, 0.1958844673505989, 0.09817479430199844, 0.04920395356814509, 0.02466039337234338])
Id16m  = np.array([0.000146803187850, 0.000151239476723, 0.0001535113826055, 0.00015509268635904, 0.00015526878069223, 0.00015556468073033])
sId16  = np.array([0.037561002693617, 0.046333193476308, 0.0353442616907872, 0.00357720125183330, 0.04892161400277213, 0.04832090850818879])
maId16 = np.array([0.065315569054276, 0.085572379349119, 0.0638077129765622, 0.00822789441806857, 0.09286890935945703, 0.09151983095756375])
Id16   = np.array([0.000138719134157, 0.000151228400257, 0.0001633065928455, 0.00015507128820955, 0.00016968842301269, 0.00014353797900319])


df4 = np.pi * 4 / (2 * N)
df8 = np.pi * 8 / (2 * N)
df16 = np.pi * 16 / (2 * N)

ax1.loglog(d4, g4, color="tab:blue", linestyle='-',marker='*', label=r"\textbf{$\mathcal{W} = 4 \omega_0$}", lw=1.2)
ax1.loglog(d8, g8, color="tab:green",  linestyle='-',marker='*', label=r"\textbf{$\mathcal{W} = 8 \omega_0$}")
ax1.loglog(d16, g16, color="tab:red",  linestyle='-',marker='*', label=r"\textbf{$\mathcal{W} = 16 \omega_0$}")

# ax1.loglog(df3, np.power(10, g3), color="tab:green", linestyle='-',marker='+', label=r"\textbf{$\mathcal{W} = 16 \omega_0$}")
ax1.loglog([5e-3, 1e1], [3.5 * 5e-3, 3.5 * 1e1], linestyle=(0, (4, 1)), color='black', lw=1)
ax1.loglog([5e-3, 1e1], [4.5 * 5e-3, 4.5 * 1e1], linestyle=(0, (1, 1)), color='black', lw=0.5)
ax1.loglog([5e-3, 1e1], [2.5 * 5e-3, 2.5 * 1e1], linestyle=(0, (1, 1)), color='black', lw=0.5)


# A = 3.3 \pm 0.5
#
ax2.loglog(d4, abs(Ig4m - IML4) / IML4, color="tab:blue", linestyle='-', marker='o', lw=1.2)
ax2.loglog(d4, abs(Id4 -  IML4) / IML4, color="tab:blue", linestyle='dashdot', marker='x', lw=1.0)
ax2.loglog(d4, abs(Id4m - IML4) / IML4, color="tab:blue", linestyle='--', marker='+', lw=1.2)
#
#ax2.loglog(d8, abs(Ig8m - IML8) / IML8, color="tab:red", linestyle='-', marker='o', lw=1.2)
#ax2.loglog(d8, abs(Id8 -  IML8) / IML8, color="tab:red", linestyle='-', marker='x', lw=1.0)
#ax2.loglog(d8, abs(Id8m - IML8) / IML8, color="tab:red", linestyle='-', marker='+', lw=1.2)
#
ax2.loglog(d16, abs(Ig16m - IML16) / IML16, color="tab:red", linestyle='-', marker='o', lw=1.2)
ax2.loglog(d16, abs(Id16 -  IML16) / IML16, color="tab:red", linestyle='dashdot', marker='x', lw=1.0)
ax2.loglog(d16, abs(Id16m - IML16) / IML16, color="tab:red", linestyle='--', marker='+', lw=1.2)

ax1.set_xlim(5e-3, 1)
ax1.set_ylim(1e-2, 3)
ax1.set_xticklabels([])
ax1.set_ylabel(r"$\check{\gamma} \ (\omega_0)$")

ax2.set_xlim(5e-3, 1)
ax2.set_ylim(1e-3, 2e-1)

ax2.set_xlabel(r"$\Delta_F \ (\omega_0)$")
ax2.set_ylabel(r"$|I_{\diamondsuit}-I_{\diamondsuit}^{M}|/|I_{\diamondsuit}^{M}|$")

ax2.loglog([1e-4, 1e-4], [1e-4, 1e-4], marker='o', linestyle='-', color = 'black', label = r"$I_\diamondsuit(\check{\gamma})$")
ax2.loglog([1e-4, 1e-4], [1e-4, 1e-4], marker='x', linestyle='dashdot', color = 'black', label = r"$I_\diamondsuit(\gamma^{\Hsquare})$")
ax2.loglog([1e-4, 1e-4], [1e-4, 1e-4], marker='+', linestyle='--', color = 'black', label = r"$\overline{I_\diamondsuit}(\gamma^\Hsquare)$")

# Create the legend
ax1.legend(loc=(0.58, 0.02),
           ncol=1,
           labelspacing=0.2,
           handlelength=1.4,
           handletextpad=0.8,
           frameon=False,
           )

ax2.legend(loc=(0.02, 0.02),
           ncol=3,
           labelspacing=0.1,
           handlelength=1.6,
           columnspacing=1.5,
           handletextpad=0.0,
           frameon=False,
           )

ax1.annotate('(a)', (-0.26, 0.89), xycoords='axes fraction', annotation_clip=False)
ax2.annotate('(b)', (-0.26, 0.98), xycoords='axes fraction', annotation_clip=False)
plt.savefig(path("Fig5.pdf"), bbox_inches='tight')