"""Literal Python conversion of fig8_20.m."""

from datetime import date as calendar_date
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np


def bbspline(x):
    # Calculate the value of a cubic B-spline at point x
    x = abs(x)
    if x > 2:
        y = 0
    else:
        if x > 1:
            y = (2 - x) ** 3 / 6
        else:
            y = 2 / 3 - x ** 2 * (1 - x / 2)
    return y


def fig8_20():
    #  DATE  TMAX  TMIN  TOBS
    #  temperature values are in Celsius degrees to tenths
    #  temp=-9999 means not recorded (need to skip these)
    #  data source: http://www.ncdc.noaa.gov
    #  date(1)=1 is 2014/01/01  and  date(end)=2015/12/16

    tempss = np.loadtxt(Path(__file__).with_name("tempss.txt"))
    AA = tempss[np.lexsort(tuple(tempss[:, j] for j in range(tempss.shape[1] - 1, -1, -1)))]
    n1, n2 = AA.shape
    n3 = 1; n4 = 1

    date_values = []
    temp_values = []
    ic = 0
    for i in range(n1):
        if AA[i, 1] != -9999:
            date_integer = int(round(AA[i, 0]))
            year = date_integer // 10000
            month = (date_integer // 100) % 100
            day = date_integer % 100
            ic = ic + 1
            datez = calendar_date(year, month, day).toordinal()
            if ic == 1:
                date0 = datez
            date_values.append(datez - date0 + 1)
            temp_values.append(AA[i, 1] * 0.1)
    date_ = np.asarray(date_values, dtype=float)
    temp = np.asarray(temp_values, dtype=float)
    nd = len(date_)
    print(f"nd = {nd}")

    ### spline fit
    m = 5
    ptsm = np.linspace(np.min(date_), np.max(date_), m)
    dp = ptsm[1] - ptsm[0]
    pts = np.concatenate(([ptsm[0] - dp], ptsm, [ptsm[m - 1] + dp]))
    B = np.zeros((m + 2, nd))
    for id_ in range(nd):
        for im in range(m + 2):
            x = (date_[id_] - pts[im]) / dp
            B[im, id_] = bbspline(x)
    C = B @ B.T
    print(f"cond(C,inf) = {np.linalg.cond(C, np.inf)}")
    d = B @ temp
    a = np.linalg.solve(C, d)
    # evaluate interpolation function s(x)
    nx = 200
    dat = np.linspace(np.min(date_), np.max(date_), nx)
    s = np.zeros(nx)
    for ix in range(nx):
        sum_ = 0
        for k in range(m + 2):
            xx = (dat[ix] - pts[k]) / dp
            sum_ = sum_ + a[k] * bbspline(xx)
        s[ix] = sum_

    plt.close("all")
    # get(gcf)
    # MATLAB Position [25 1115 658 230], scaled to fit the plot window.
    fig, ax = plt.subplots(figsize=(6.58, 3.0), layout="constrained")
    co = np.array([
        [0, 0, 1],
        [0, 0.5, 0],
        [1, 0, 0],
        [0, 0.75, 0.75],
        [0.75, 0, 0.75],
        [0.75, 0.75, 0],
        [0.25, 0.25, 0.25],
    ])
    ax.set_prop_cycle(color=co)

    ### plot entire data set
    # ax.plot(date_,temp,'or',markersize=8,linewidth=2)

    ### plot nx (randomly choosen) data points
    nx = 50
    r = np.random.randint(0, nd, size=nx)
    d = np.zeros(nx)
    t = np.zeros(nx)
    for i in range(nx):
        d[i] = date_[r[i]]
        t[i] = temp[r[i]]

    ax.plot(d, t, "or", markersize=8, linewidth=2)
    ax.plot(dat, s, "b", linewidth=2)
    ax.grid(True)
    ax.set_frame_on(True)
    ax.set_xlabel("Date", fontsize=12, fontweight="bold")
    ax.set_ylabel("Temperature", fontsize=12, fontweight="bold")
    tmax = calendar_date(2015, 12, 31).toordinal() - date0 + 1
    print(f"tmax = {tmax}")
    t2 = calendar_date(2014, 7, 1).toordinal() - date0 + 1
    t3 = calendar_date(2015, 1, 1).toordinal() - date0 + 1
    t4 = calendar_date(2015, 7, 1).toordinal() - date0 + 1
    ax.axis([0, tmax, -20, 40])
    ax.set_xticks([1, t2, t3, t4, tmax])
    ax.set_xticklabels(["Jan", "July", "Jan", "July", "Jan"])
    ax.tick_params(labelsize=12)
    for label in ax.get_xticklabels() + ax.get_yticklabels():
        label.set_fontweight("bold")

    ### generate data for exercise in text
    nxx = 15
    dd = np.zeros(nxx)
    tt = np.zeros(nxx)
    dd[0] = date_[0]
    tt[0] = temp[0]
    dd[nxx - 1] = date_[nd - 1]
    tt[nxx - 1] = temp[nd - 1]
    print(nd)
    nn = round(nd / (nxx - 1))
    print(nn)
    i1 = 0
    for i in range(nn - 1, nd, nn):
        i1 = i1 + 1
        iz = i + np.random.randint(1, 17) - 8
        dd[i1] = date_[iz]
        tt[i1] = temp[iz]
    print(dd)
    print(tt)

    # fig.savefig('/Users/mark/Desktop/temp.eps')

    # fig2, ax2 = plt.subplots(figsize=(6.58,3.0), layout='constrained')
    # ax2.axis([0,tmax,-20,40])
    # ax2.plot(dd,tt,'ro',linewidth=2)

    plt.show()


if __name__ == "__main__":
    fig8_20()
