import numpy as np
import scipy.fftpack
import zlib

# Generate a sample sequence of elapsed time (in nanoseconds)
#sequence = np.random.randint(1, 1000, size=1000)  # Example sequence
#sequence = [0,0,0,0,50000,35000,40000,20000,35000,300,400,200,400,500,300,200,0,0,0]
sequence=[5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          5000,4000,5500,4500,3000,4000,
          1000000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
          2500,3000,2500,4000,
         ]


# Apply the DCT
dct_coefficients = scipy.fftpack.dct(sequence, norm='ortho')

# Zero out mid-frequency components
low_cutoff = int(0.3 * len(dct_coefficients))
high_cutoff = int(0.7 * len(dct_coefficients))
quantized_coefficients = dct_coefficients.copy()
#quantized_coefficients[low_cutoff:high_cutoff] = 0
quantized_coefficients[5:-5] =0
#print(quantized_coefficients)

# Compress the quantized coefficients using zlib
compressed_data = zlib.compress(quantized_coefficients.astype(np.float32).tobytes())
print('compressed len',len(compressed_data))

# Decompress the data
decompressed_data = np.frombuffer(zlib.decompress(compressed_data), dtype=np.float32)

# Apply the inverse DCT
reconstructed_sequence = scipy.fftpack.idct(decompressed_data, norm='ortho')

# Post-process the reconstructed sequence
# Ensure non-negative values
reconstructed_sequence[reconstructed_sequence < 0] = 0

# Scale the reconstructed sequence to ensure the sum matches the original
scaling_factor = np.sum(sequence) / np.sum(reconstructed_sequence)
reconstructed_sequence *= scaling_factor

# Convert to integer values, if necessary
reconstructed_sequence = np.round(reconstructed_sequence).astype(int)

# Ensure the sum of gaps is exactly the same as the original
adjustment = np.sum(sequence) - np.sum(reconstructed_sequence)
reconstructed_sequence[-1] += adjustment

print("Original sequence sum:", np.sum(sequence))
print("Reconstructed sequence sum:", np.sum(reconstructed_sequence))

# Plot the original and reconstructed sequences for comparison
import matplotlib.pyplot as plt

plt.figure(figsize=(12, 6))
plt.subplot(2, 1, 1)
plt.plot(sequence, label='Original Sequence')
plt.legend()
plt.subplot(2, 1, 2)
plt.plot(reconstructed_sequence, label='Reconstructed Sequence')
plt.legend()
plt.show()