"""Hintereisferner with centerlines: a run from the v5 `centerlines` gdir.

To check how much of the length misfit of lam8_jjas_dyn (elevation bands)
comes from the flowline representation. The v5 centerlines gdirs are the
"default" experiment (trapezoid_lambdas 2, consensus, prcp_fac 2.5, melt_f
5) with centerlines, statically calibrated, and were built with the old
jjas_hists_v20261001 forcing. So here:
- the climate is re-processed from jjas_hists_v20261005 ens183, the forcing
  of lam8_jjas_dyn. The two are identical from 1900, so the static
  calibration (Hugonnet 2000-2020) stays valid: this is asserted.
- the run is the round 4 one: from the RGI-date flowlines in year 1, with
  evolution_model FluxBased (as the centerlines runs of round 4).
- no dynamic melt_f calibration.

Its own working dir (_wd_centerlines), nothing is shared with the elevation
band gdir of validate_hef.py. Output: hef_centerlines_run.csv (volume, area,
length per year), read by validate_hef.py.

Needs lmr_forcing on the path:
PYTHONPATH=../../lmr_gdirs_v5_202610 python run_hef_centerlines.py
"""
import os

import numpy as np
import xarray as xr

import oggm.cfg as cfg
from oggm import utils, workflow, tasks

import lmr_forcing

RGI_ID = 'RGI60-11.00897'
GDIR_EXP = 'centerlines'
HIST_DIR = 'jjas_hists_v20261005'
MEMBER = 183
SUFFIX = f'_{HIST_DIR}_ens{MEMBER:03d}'

HERE = os.path.dirname(os.path.abspath(__file__))
V5_EXP_DIR = os.path.join(HERE, '..', '..', 'lmr_gdirs_v5_202610', 'exps',
                          GDIR_EXP)
BASE_URL = 'https://cluster.klima.uni-bremen.de/~fmaussion/runs/lmr_2025/'
GDIR_URL = BASE_URL + f'lmr_gdirs_v5_202610/exps/{GDIR_EXP}/'
LMR_URL = BASE_URL + 'data/unstacked/lmr_seas_fullsplit_oldBE'

cfg.initialize(file=os.path.join(V5_EXP_DIR, 'params.cfg'),
               logging_level='WARNING')
cfg.PARAMS['border'] = 160
cfg.PATHS['working_dir'] = utils.mkdir(os.path.join(HERE, '_wd_centerlines'))
# As the round 4 runs
cfg.PARAMS['store_diagnostic_variables'] = ['volume', 'volume_bsl', 'area',
                                            'area_min_h', 'length']
cfg.PARAMS['min_ice_thick_for_area'] = 2.
cfg.PARAMS['min_ice_thick_for_length'] = 2.
cfg.PARAMS['evolution_model'] = 'FluxBased'

gdir = workflow.init_glacier_directories([RGI_ID], from_prepro_level=3,
                                         prepro_base_url=GDIR_URL,
                                         prepro_rgi_version='62',
                                         prepro_border=160)[0]
print(f'{RGI_ID} centerlines: melt_f {gdir.settings["melt_f"]}, prcp_fac '
      f'{gdir.settings["prcp_fac"]:.2f}, temp_bias '
      f'{gdir.settings["temp_bias"]:.2f}, '
      f'{len(gdir.read_store("model_flowlines"))} flowlines')

lmr_forcing.process_mira_data(gdir, ens_member=MEMBER, server=LMR_URL,
                              temp_var='tp_monthly',
                              temp_server=LMR_URL + '/' + HIST_DIR,
                              output_filesuffix=SUFFIX)

# The new climate must be the gdir's own from 1900 (the calibration period)
tc = xr.coders.CFDatetimeCoder(use_cftime=True)
with xr.open_dataset(gdir.get_filepath('climate_historical'),
                     decode_times=tc) as old, \
        xr.open_dataset(gdir.get_filepath('climate_historical',
                                          filesuffix=SUFFIX),
                        decode_times=tc) as new:
    for v in ['temp', 'prcp']:
        d = np.abs(old[v].values - new[v].values)
        yr = old['time.year'].values
        print(f'{v}: max |new - gdir| before 1850 {d[yr < 1850].max():.4f}, '
              f'1850-1899 {d[(yr >= 1850) & (yr < 1900)].max():.4f}, '
              f'from 1900 {d[yr >= 1900].max():.6f}')
        assert d[yr >= 1900].max() < 1e-4

tasks.run_from_climate_data(gdir, ys=1,
                            climate_filename='climate_historical',
                            climate_input_filesuffix=SUFFIX,
                            output_filesuffix='_hef_centerlines')
with xr.open_dataset(gdir.get_filepath('model_diagnostics',
                                       filesuffix='_hef_centerlines')) as ds:
    df = ds[['volume_m3', 'area_m2', 'length_m']].to_dataframe()
df = df[['volume_m3', 'area_m2', 'length_m']]
df.columns = ['volume', 'area', 'length']
df.index = df.index.astype(int)
df.index.name = 'year'
df.to_csv(os.path.join(HERE, 'hef_centerlines_run.csv'))
print(df.loc[[1, 1700, 1850, 1900, 1950, 2003, 2020, df.index[-1]]])