#include <algorithm>
#include <cpgplot.h>
#include <iostream>
#include <random>
#include <set>
#include <string>
#include <vector>

using namespace std;

struct Color {
    struct ColRs {
        float lower, upper;
        int clcount;
        float Rl, Rr;
        float Gl, Gr;
        float Bl, Br;
        bool operator< (ColRs other) const {
            return lower < other.lower;
        }
    };
    int start;
    vector<pair<float, int>> table;
    static float dist(int it, int nall, float start, float end) {
        return start + (end - start) * it / (nall - 1);
    }
    int getidx(float val) {
        auto it = lower_bound(table.begin(), table.end(), pair{val, 0})
         - table.begin();
        if (it && (it == table.size() ||
                val - table[it - 1].first <= table[it].first - val))
            return start + it - 1;
        return start + it;
    }
    Color(const set<ColRs>& ctable, int start_id)
        : start(start_id) {
            for (auto i : ctable)
            {
                for (int j = 0; j < i.clcount; j++)
                {
                    cpgscr(start_id,
                        dist(j, i.clcount, i.Rl, i.Rr) / 255.,
                        dist(j, i.clcount, i.Gl, i.Gr) / 255.,
                        dist(j, i.clcount, i.Bl, i.Br) / 255.
                        );
                    table.emplace_back(dist(j, i.clcount, i.lower, i.upper), start_id);
                    ++start_id;
                }
            }
        cpgscir(start, start_id - 1);
    }
};

struct PGgrid{
    void print(const int* grid, int height, int width) {
        cpgpixl(grid, height, width, 1, height, 1, width,
                0, 1, 0, 1);
    }
    PGgrid(std::string device = "/XWINDOW") {
        if (cpgopen(device.c_str()) < 0)
            throw std::logic_error("PGgrid(" + device + ")");
    }
    ~PGgrid(){ cpgend(); }
};

int main() {
    PGgrid grid;
    int N = 100, M = 100;
    float arr[N * M];
    Color rammb({
        {-105, -90, 10, 56, 255, 56, 255, 56, 255}, // 白 T=-105...-90,R=G=B255->56
        {-90, -80, 10, 90, 255, 90, 255, 77, 84}, // 金 T=-90...-80,R=G255->90,B84->77
        {-80, -70, 10, 226, 84, 28, 49, 8, 34}, // 红 T=-80...-70,R84->226,G28->49,B8->34
        {-70, -60, 10, 113, 41, 242, 98, 73, 24}, // 绿 T=-70...-60,R41->113,G98->242,B24->73
        {-60, -50, 10, 0, 0, 0, 0, 236, 96}, // 蓝 T=-60...-50,R0,G0,B96->236
        {-50, -30, 10, 88, 196, 90, 255, 90, 255}, // 青 T=-50...-30,R196->88,G255->90,B255->90
        {-30, 40, 24, 255, 0, 255, 0, 255, 0} // 灰度（T=-30...40,R=G=B0->255）
    }, 16);
    cpgscf(2);
    cpglab("(D_1)", "(D_2)", "Test I");
    for (int i = 0; i < N; i++)
        for (int j = 0; j < M; j++)
            arr[i * M + j] = sqrt(pow(i - N / 2, 2) + pow(j - M / 2, 2))
        / sqrt(pow(N / 2, 2) + pow(M / 2, 2)) * 140 - 100;
    int pst[N * M];
    for (int i = 0; i < N * M; i++)
        pst[i] = rammb.getidx(arr[i]);
    grid.print(pst, N, M);
    cpgopen("/XWIN");
    float fixv[] = {-80, -70, -60, -50, -40, -30, 10}, tr[] = {0, 1.0f / N, 0, 0, 0, 1.0f / M};
    cpgcont(arr, N, M, 1, N, 1, M, fixv, 7, tr);
    return 0;
}

#include <cpgplot.h>
#include <iostream>
#include <random>

void init_BD() {
    int col_id = 16;
    auto newcol = [&] (float GrD) {
        cpgscr(col_id++, GrD / 255, GrD / 255, GrD / 255);
    };
    newcol(85); // -84
    newcol(135); // -81
    newcol(135); // -78
    newcol(255); // -75
    newcol(255); // -72
    newcol(0); // -69
    newcol(0); // -66
    newcol(160); // -63
    newcol(160); // -60
    newcol(160); // -57
    newcol(160); // -54
    newcol(110); // -51
    newcol(110); // -48
    newcol(110); // -45
    newcol(110); // -42
    newcol(60); // -39
    newcol(60); // -36
    newcol(60); // -33
    newcol(198.750000); // -30
    newcol(192.000000); // -27
    newcol(185.250000); // -24
    newcol(178.500000); // -21
    newcol(171.750000); // -18
    newcol(165.000000); // -15
    newcol(158.250000); // -12
    newcol(151.500000); // -9
    newcol(144.750000); // -6
    newcol(138.000000); // -3
    newcol(131.250000); // 0
    newcol(124.500000); // 3
    newcol(117.750000); // 6
    newcol(111.000000); // 9
    newcol(213.428571); // 12
    newcol(177.857143); // 15
    newcol(142.285714); // 18
    newcol(106.714286); // 21
    newcol(71.142857); // 24
    newcol(35.571429); // 27
    newcol(0.000000); // 30
    cpgscir(16, col_id - 1);
}

int main() {
    using namespace std;
    cpgopen("/XWINDOW");
    init_BD();
    random_device seed;
    mt19937 rnd(seed());
    int H = 1000, W = 1600;
    float arr[H * W];
    uniform_real_distribution<float> dist(-90, 30);
    for (int i = 0; i < H; i++)
        for (int j = 0; j < W; j++)
            arr[i * W + j] = dist(rnd);
    float trans[6]{-0.5, 0, 1, -0.5, 1, 0};
    cpgenv(0, W, 0, H, 1, 1);
    cpglab("LON", "LAT", "IR-BD");
    cpgimag(arr, H, W, 1, H, 1, W, -84, 30, trans);
    cpgwedg("RI", 1, 3, -84, 30, "\20C");
    cpgend();
    return 0;
}