#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Fri Jun  5 11:53:53 2026

@author: leonhardtm
"""

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from pathlib import Path
import qutip as qt

##### (a) Mode Occupations

M=pd.read_csv(Path(__file__).resolve().parent.parent/"DATA"/"FIG4"/'ModeOccupations.csv')
M_data=M.to_numpy()

plt.figure()
plt.plot(M_data[:,0],M_data[:,1],M_data[:,0],M_data[:,2])
legend=plt.legend(['1','2'],title='Mode')
plt.xlabel('t')
plt.ylabel(r'$\langle\hat{a}^\dagger_k\hat{a}_k\rangle/N$')
plt.xlim([0,1500])


##### (b) Indivdual Cluster Trajectories and Periodic Orbits

for i in range(3):
    T=pd.read_csv(Path(__file__).resolve().parent.parent/"DATA"/"FIG4"/'Trajectory_{}.csv'.format(i+1))
    T_data=T.to_numpy()
    O=pd.read_csv(Path(__file__).resolve().parent.parent/"DATA"/"FIG4"/'Orbit_{}.csv'.format(i+1))
    O_data=O.to_numpy()
    
    b=qt.Bloch()
    b.view=[-60,30]
    b.add_points([2*T_data[:,1],2*T_data[:,2],2*T_data[:,3]],meth='l')
    b.add_points([2*O_data[:,0],2*O_data[:,1],2*O_data[:,2]],meth='l')
    b.show()
    b.axes.set_title('Cluster {}'.format(i+1))