
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
from scipy.stats import linregress
from scipy.interpolate import splrep, splev
import requests
from io import StringIO

# 1. Data Acquisition
# --- NEW URL for Antarctic Mass Data ---
url = 'https://archive.podaac.earthdata.nasa.gov/podaac-ops-cumulus-protected/ANTARCTICA_MASS_TELLUS_MASCON_CRI_TIME_SERIES_RL06.3_V4/antarctica_mass_200204_202502.txt'

# 2. Fetch the data using requests with .netrc authentication
try:
    response = requests.get(url)
    response.raise_for_status()
    data_content = response.text
except requests.exceptions.RequestException as e:
    print(f"Error fetching data: {e}")
    print("Please ensure you have set up your NASA Earthdata Login and .netrc file correctly.")
    exit()

# Use io.StringIO to treat the string content as a file
data_io = StringIO(data_content)

# Define column names for Antarctic Mass data based on file header
# --- UPDATED Column Names (3 columns, similar to Greenland) ---
column_names = [
    'Decimal_Year',                  # Col 1: TIME (year.decimal)
    'Antarctica_Mass_Gt',            # Col 2: Antarctica mass (Gigatonnes) - Our Y-axis
    'Antarctica_Mass_Uncertainty_Gt' # Col 3: Antarctica mass 1-sigma uncertainty (Gigatonnes)
]

# 3. Read the data into a Pandas DataFrame
df = pd.read_csv(
    data_io,
    sep=r'\s+',
    engine='python',
    comment='#',
    # --- UPDATED skiprows for Antarctic Mass data (expected to be 34) ---
    skiprows=34, # Expected based on similar files (Greenland), adjust if needed
    names=column_names,
)
df


# 4. Prepare data for plotting
# --- UPDATED to use Antarctica_Mass_Gt for Y-data ---
sea_level_column = 'Antarctica_Mass_Gt'
y_axis_label = 'Antarctica Mass Anomaly (Gigatonnes)' # Updated Y-axis label

df.dropna(subset=[sea_level_column], inplace=True) # Drop NaNs if any are present

# Extract data for plotting and spline fitting
x_data = df['Decimal_Year'].to_numpy()
y_data = df[sea_level_column].to_numpy()

# --- Cubic Spline Fitting ---
spl_representation = splrep(x_data, y_data, k=3, s=1000000) # Retaining user's preferred high smoothing
x_spline = np.linspace(x_data.min(), x_data.max(), 500)
y_spline = splev(x_spline, spl_representation)
# ------------------------------------------

# --- Linear Trend Calculation (unchanged logic) ---
original_slope, original_intercept, _, _, _ = linregress(x_data, y_data)
trend_line = original_slope * x_data + original_intercept
# ----------------------------------------------------

# 5. Plot the time series
plt.figure(figsize=(14, 7))

# Plot the original data
plt.plot(x_data, y_data, color='#4FC3F7', linestyle='-', label='Antarctica Mass Anomaly', alpha=0.8)

# Plot the main linear trend line
plt.plot(x_data, trend_line, color='red', linestyle='--', label=f'Linear Trend: {original_slope:.2f} Gt/year')

# Plot the Cubic Spline line
plt.plot(x_spline, y_spline, color='#673AB7', linestyle='-', linewidth=2, label='Cubic Spline Fit')

# 7. Customize the plot
# --- UPDATED Plot Title and Figtext Source ---
plt.title('Antarctica Ice Sheet Mass Change (NASA GRACE/GRACE-FO Satellite Data)')
plt.xlabel('Year')
plt.ylabel(y_axis_label)
plt.grid(True, linestyle=':', alpha=0.7)
plt.legend()
plt.tight_layout()

plt.figtext(0.8, 0.01, "Source: NASA PO.DAAC (ANTARCTICA_MASS_TELLUS_MASCON_CRI_TIME_SERIES_RL06.3_V4)\nColumn 2: Antarctica mass (Gigatonnes)", ha="center", fontsize=9, color="gray")

plt.show()

print(f"\nOriginal Linear Trend Slope (Antarctica Mass): {original_slope:.2f} Gt/year")
print(f"R-squared value for linear fit (Antarctica Mass): {linregress(x_data, y_data).rvalue**2:.4f}")

Original Linear Trend Slope (Antarctica Mass): -133.81 Gt/year
R-squared value for linear fit (Antarctica Mass): 0.9512

	Decimal_Year	Antarctica_Mass_Gt	Antarctica_Mass_Uncertainty_Gt
0	2002.71	71.70	97.10
1	2002.79	85.80	61.76
2	2002.87	-48.73	58.63
3	2002.96	-2.54	59.14
4	2003.04	13.34	64.66
...	...	...	...
234	2024.79	-2370.20	47.19
235	2024.87	-2597.23	58.88
236	2024.96	-2701.84	72.07
237	2025.04	-2636.40	85.84
238	2025.12	-2655.00	102.01