"""
This plot presents a 2-D histogram of reported cloud ceilings as inferred by
the METAR cloud amount reports of at least BKN (broken) or OVC (overcast)
from the automated sensor.  Please note that the yaxis
uses an irregular spacing.
"""

from datetime import date

import matplotlib.colors as mpcolors
import numpy as np
import numpy.ma as ma
import pandas as pd
from pyiem.database import get_sqlalchemy_conn, sql_helper
from pyiem.exceptions import NoDataFound
from pyiem.plot import figure_axes, get_cmap


def get_description():
    """Return a dict describing how to call this plotter"""
    desc = {"description": __doc__, "data": True, "cache": 86400}
    desc["arguments"] = [
        dict(
            type="zstation",
            name="zstation",
            default="AMW",
            network="IA_ASOS",
            label="Select Station:",
        ),
        dict(type="cmap", name="cmap", default="jet", label="Color Ramp:"),
    ]
    return desc


def plotter(ctx: dict):
    """Go"""
    station = ctx["zstation"]
    with get_sqlalchemy_conn("asos") as conn:
        df = pd.read_sql(
            sql_helper("""
            select extract(doy from valid) as doy,
            CASE
                WHEN skyc1 in ('BKN', 'OVC') THEN skyl1
                WHEN skyc2 in ('BKN', 'OVC') THEN skyl2
                WHEN skyc3 in ('BKN', 'OVC') THEN skyl3
                WHEN skyc4 in ('BKN', 'OVC') THEN skyl4
            END as sky
            from alldata
            WHERE station = :station
            and valid > '1973-01-01' and report_type = 3
            and (skyc1 in ('BKN', 'OVC') or skyc2 in ('BKN', 'OVC')
                 or skyc3 in ('BKN', 'OVC') or skyc4 in ('BKN', 'OVC'))
        """),
            conn,
            params={"station": station},
            index_col=None,
        )
    if df.empty:
        raise NoDataFound("Error, no results returned!")

    w = np.arange(1, 366, 7)
    z = np.append(
        np.arange(100, 5001, 100),
        np.append(
            np.arange(5500, 10001, 500),
            np.arange(11000, 31001, 1000),
        ),
    )

    H, xedges, yedges = np.histogram2d(
        df["sky"].to_numpy(), df["doy"].to_numpy(), bins=(z, w)
    )
    rows = []
    for i, x in enumerate(xedges[:-1]):
        for j, y in enumerate(yedges[:-1]):
            rows.append(dict(ceiling=x, doy=y, count=H[i, j]))
    resdf = pd.DataFrame(rows)

    H = ma.array(H)
    H.mask = np.where(H < 1, True, False)

    ab = ctx["_nt"].sts[station]["archive_begin"]
    if ab is None:
        raise NoDataFound("Unknown station metadata.")
    syear = max([1973, ab.year])

    (fig, ax) = figure_axes(
        title=f"({syear}-{date.today().year}) {ctx['_sname']} Ceiling",
        subtitle="Frequency of first cloud level with BKN or OVC reported",
        apctx=ctx,
    )

    bounds = np.arange(0, 1.2, 0.1)
    bounds = np.concatenate((bounds, np.arange(1.2, 2.2, 0.2)))
    cmap = get_cmap(ctx["cmap"])
    cmap.set_under("#F9CCCC")
    # If the number of bounds is greater than the cmap colors, we need to
    # adjust
    if len(bounds) > cmap.N:
        bounds = np.linspace(0, 2, cmap.N)
    norm = mpcolors.BoundaryNorm(bounds, cmap.N)

    years = date.today().year - syear + 1.0
    c = ax.imshow(
        H / years, aspect="auto", interpolation="nearest", norm=norm, cmap=cmap
    )
    ax.set_ylim(-0.5, len(z) - 0.5)
    idx = [0, 4, 9, 19, 29, 39, 49, 54, 59, 64, 69, 74, 79]
    ax.set_yticks(idx)
    ax.set_yticklabels(z[idx])
    ax.set_ylabel("Overcast/Broken Level [ft AGL], irregular scale")
    ax.set_xlabel("Week of the Year")
    ax.set_xticks(np.arange(1, 55, 7))
    ax.set_xticklabels(
        "Jan 1,Feb 19,Apr 8,May 27,Jul 15,Sep 2,Oct 21,Dec 9".split(",")
    )
    b = fig.colorbar(c)
    b.set_label("Hourly Obs per week per year")
    return fig, resdf
