Actual source code: cmap.c

      1: #include <petscsys.h>
      2: #include <petscdraw.h>

      4: /*
      5:     Set up a color map, using uniform separation in hue space.
      6:     Map entries are Red, Green, Blue.
      7:     Values are "gamma" corrected.
      8:  */

     10: /*
     11:    Gamma is a monitor dependent value.  The value here is an
     12:    approximate that gives somewhat better results than Gamma = 1.
     13:  */
     14: static PetscReal Gamma = 2.0;

     16: PetscErrorCode PetscDrawUtilitySetGamma(PetscReal g)
     17: {
     18:   PetscFunctionBegin;
     19:   Gamma = g;
     20:   PetscFunctionReturn(PETSC_SUCCESS);
     21: }

     23: static inline double PetscHlsHelper(double m1, double m2, double h)
     24: {
     25:   while (h > 1.0) h -= 1.0;
     26:   while (h < 0.0) h += 1.0;
     27:   if (h < 1 / 6.0) return m1 + (m2 - m1) * h * 6;
     28:   if (h < 1 / 2.0) return m2;
     29:   if (h < 2 / 3.0) return m1 + (m2 - m1) * (2 / 3.0 - h) * 6;
     30:   return m1;
     31: }

     33: static inline void PetscHlsToRgb(double h, double l, double s, double *r, double *g, double *b)
     34: {
     35:   if (s > 0.0) {
     36:     double m2 = l <= 0.5 ? l * (1.0 + s) : l + s - (l * s);
     37:     double m1 = 2 * l - m2;
     38:     *r        = PetscHlsHelper(m1, m2, h + 1 / 3.);
     39:     *g        = PetscHlsHelper(m1, m2, h);
     40:     *b        = PetscHlsHelper(m1, m2, h - 1 / 3.);
     41:   } else {
     42:     /* ignore hue */
     43:     *r = *g = *b = l;
     44:   }
     45: }

     47: static inline void PetscGammaCorrect(double *r, double *g, double *b)
     48: {
     49:   PetscReal igamma = 1 / Gamma;
     50:   *r               = (double)PetscPowReal((PetscReal)*r, igamma);
     51:   *g               = (double)PetscPowReal((PetscReal)*g, igamma);
     52:   *b               = (double)PetscPowReal((PetscReal)*b, igamma);
     53: }

     55: static PetscErrorCode PetscDrawCmap_Hue(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
     56: {
     57:   int    i;
     58:   double maxhue = 212.0 / 360, lightness = 0.5, saturation = 1.0;

     60:   PetscFunctionBegin;
     61:   for (i = 0; i < mapsize; i++) {
     62:     double hue = maxhue * (double)i / (mapsize - 1), r, g, b;
     63:     PetscHlsToRgb(hue, lightness, saturation, &r, &g, &b);
     64:     PetscGammaCorrect(&r, &g, &b);
     65:     R[i] = (unsigned char)(255 * PetscMin(r, 1.0));
     66:     G[i] = (unsigned char)(255 * PetscMin(g, 1.0));
     67:     B[i] = (unsigned char)(255 * PetscMin(b, 1.0));
     68:   }
     69:   PetscFunctionReturn(PETSC_SUCCESS);
     70: }

     72: static PetscErrorCode PetscDrawCmap_Gray(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
     73: {
     74:   PetscFunctionBegin;
     75:   for (int i = 0; i < mapsize; i++) R[i] = G[i] = B[i] = (unsigned char)((255.0 * i) / (mapsize - 1));
     76:   PetscFunctionReturn(PETSC_SUCCESS);
     77: }

     79: static PetscErrorCode PetscDrawCmap_Jet(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
     80: {
     81:   int          i;
     82:   const double knots[] = {0, 1 / 8., 3 / 8., 5 / 8., 7 / 8., 1};

     84:   PetscFunctionBegin;
     85:   for (i = 0; i < mapsize; i++) {
     86:     double u = (double)i / (mapsize - 1);
     87:     double m, r = 0, g = 0, b = 0;
     88:     int    k = 0;
     89:     while (k < 4 && u > knots[k + 1]) k++;
     90:     m = (u - knots[k]) / (knots[k + 1] - knots[k]);
     91:     switch (k) {
     92:     case 0:
     93:       r = 0;
     94:       g = 0;
     95:       b = (m + 1) / 2;
     96:       break;
     97:     case 1:
     98:       r = 0;
     99:       g = m;
    100:       b = 1;
    101:       break;
    102:     case 2:
    103:       r = m;
    104:       g = 1;
    105:       b = 1 - m;
    106:       break;
    107:     case 3:
    108:       r = 1;
    109:       g = 1 - m;
    110:       b = 0;
    111:       break;
    112:     case 4:
    113:       r = 1 - m / 2;
    114:       g = 0;
    115:       b = 0;
    116:       break;
    117:     }
    118:     R[i] = (unsigned char)(255 * PetscMin(r, 1.0));
    119:     G[i] = (unsigned char)(255 * PetscMin(g, 1.0));
    120:     B[i] = (unsigned char)(255 * PetscMin(b, 1.0));
    121:   }
    122:   PetscFunctionReturn(PETSC_SUCCESS);
    123: }

    125: static PetscErrorCode PetscDrawCmap_Hot(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
    126: {
    127:   int          i;
    128:   const double knots[] = {0, 3 / 8., 3 / 4., 1};

    130:   PetscFunctionBegin;
    131:   for (i = 0; i < mapsize; i++) {
    132:     double u = (double)i / (mapsize - 1);
    133:     double m, r = 0, g = 0, b = 0;
    134:     int    k = 0;
    135:     while (k < 2 && u > knots[k + 1]) k++;
    136:     m = (u - knots[k]) / (knots[k + 1] - knots[k]);
    137:     switch (k) {
    138:     case 0:
    139:       r = m;
    140:       g = 0;
    141:       b = 0;
    142:       break;
    143:     case 1:
    144:       r = 1;
    145:       g = m;
    146:       b = 0;
    147:       break;
    148:     case 2:
    149:       r = 1;
    150:       g = 1;
    151:       b = m;
    152:       break;
    153:     }
    154:     R[i] = (unsigned char)(255 * PetscMin(r, 1.0));
    155:     G[i] = (unsigned char)(255 * PetscMin(g, 1.0));
    156:     B[i] = (unsigned char)(255 * PetscMin(b, 1.0));
    157:   }
    158:   PetscFunctionReturn(PETSC_SUCCESS);
    159: }

    161: static PetscErrorCode PetscDrawCmap_Bone(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
    162: {
    163:   PetscFunctionBegin;
    164:   (void)PetscDrawCmap_Hot(mapsize, R, G, B);
    165:   for (int i = 0; i < mapsize; i++) {
    166:     double u = (double)i / (mapsize - 1);
    167:     double r = (7 * u + B[i] / 255.0) / 8;
    168:     double g = (7 * u + G[i] / 255.0) / 8;
    169:     double b = (7 * u + R[i] / 255.0) / 8;
    170:     R[i]     = (unsigned char)(255 * PetscMin(r, 1.0));
    171:     G[i]     = (unsigned char)(255 * PetscMin(g, 1.0));
    172:     B[i]     = (unsigned char)(255 * PetscMin(b, 1.0));
    173:   }
    174:   PetscFunctionReturn(PETSC_SUCCESS);
    175: }

    177: #include "../src/sys/classes/draw/utils/cmap/coolwarm.h"
    178: #include "../src/sys/classes/draw/utils/cmap/parula.h"
    179: #include "../src/sys/classes/draw/utils/cmap/viridis.h"
    180: #include "../src/sys/classes/draw/utils/cmap/plasma.h"
    181: #include "../src/sys/classes/draw/utils/cmap/inferno.h"
    182: #include "../src/sys/classes/draw/utils/cmap/magma.h"

    184: static struct {
    185:   const char *name;
    186:   const unsigned char (*data)[3];
    187:   PetscErrorCode (*cmap)(int, unsigned char[], unsigned char[], unsigned char[]);
    188: } PetscDrawCmapTable[] = {
    189:   {"hue",      NULL,                   PetscDrawCmap_Hue }, /* varying hue with constant lightness and saturation */
    190:   {"gray",     NULL,                   PetscDrawCmap_Gray}, /* black to white with shades of gray */
    191:   {"bone",     NULL,                   PetscDrawCmap_Bone}, /* black to white with gray-blue shades */
    192:   {"jet",      NULL,                   PetscDrawCmap_Jet }, /* rainbow-like colormap from NCSA, University of Illinois */
    193:   {"hot",      NULL,                   PetscDrawCmap_Hot }, /* black-body radiation */
    194:   {"coolwarm", PetscDrawCmap_coolwarm, NULL              }, /* ParaView default (Cool To Warm with Diverging interpolation) */
    195:   {"parula",   PetscDrawCmap_parula,   NULL              }, /* MATLAB (default since R2014b) */
    196:   {"viridis",  PetscDrawCmap_viridis,  NULL              }, /* matplotlib 1.5 (default since 2.0) */
    197:   {"plasma",   PetscDrawCmap_plasma,   NULL              }, /* matplotlib 1.5 */
    198:   {"inferno",  PetscDrawCmap_inferno,  NULL              }, /* matplotlib 1.5 */
    199:   {"magma",    PetscDrawCmap_magma,    NULL              }, /* matplotlib 1.5 */
    200: };

    202: PetscErrorCode PetscDrawUtilitySetCmap(const char colormap[], int mapsize, unsigned char R[], unsigned char G[], unsigned char B[])
    203: {
    204:   int         i, j;
    205:   const char *cmap_name_list[PETSC_STATIC_ARRAY_LENGTH(PetscDrawCmapTable)];
    206:   PetscInt    id = 0, count = PETSC_STATIC_ARRAY_LENGTH(cmap_name_list);
    207:   PetscBool   reverse = PETSC_FALSE, brighten = PETSC_FALSE;
    208:   PetscReal   beta = 0;

    210:   PetscFunctionBegin;
    211:   for (i = 0; i < count; i++) cmap_name_list[i] = PetscDrawCmapTable[i].name;
    212:   if (colormap && colormap[0]) {
    213:     PetscBool match = PETSC_FALSE;
    214:     for (id = 0; !match && id < count; id++) PetscCall(PetscStrcasecmp(colormap, cmap_name_list[id], &match));
    215:     PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Colormap '%s' not found", colormap);
    216:   }
    217:   PetscCall(PetscOptionsGetEList(NULL, NULL, "-draw_cmap", cmap_name_list, count, &id, NULL));
    218:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-draw_cmap_reverse", &reverse, NULL));
    219:   PetscCall(PetscOptionsGetReal(NULL, NULL, "-draw_cmap_brighten", &beta, &brighten));
    220:   PetscCheck(!brighten || (beta > (PetscReal)-1 && beta < (PetscReal)+1), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "brighten parameter %g must be in the range (-1,1)", (double)beta);

    222:   if (PetscDrawCmapTable[id].cmap) {
    223:     PetscCall(PetscDrawCmapTable[id].cmap(mapsize, R, G, B));
    224:   } else {
    225:     const unsigned char (*rgb)[3] = PetscDrawCmapTable[id].data;
    226:     PetscCheck(mapsize == 256 - PETSC_DRAW_BASIC_COLORS, PETSC_COMM_SELF, PETSC_ERR_SUP, "Colormap '%s' with size %d not supported", cmap_name_list[id], mapsize);
    227:     for (i = 0; i < mapsize; i++) {
    228:       R[i] = rgb[i][0];
    229:       G[i] = rgb[i][1];
    230:       B[i] = rgb[i][2];
    231:     }
    232:   }

    234:   if (reverse) {
    235:     i = 0;
    236:     j = mapsize - 1;
    237:     while (i < j) {
    238: #define SWAP(a, i, j) \
    239:   do { \
    240:     unsigned char t = a[i]; \
    241:     a[i]            = a[j]; \
    242:     a[j]            = t; \
    243:   } while (0)
    244:       SWAP(R, i, j);
    245:       SWAP(G, i, j);
    246:       SWAP(B, i, j);
    247: #undef SWAP
    248:       i++;
    249:       j--;
    250:     }
    251:   }

    253:   if (brighten) {
    254:     PetscReal gamma = (beta > 0.0) ? (1 - beta) : (1 / (1 + beta));
    255:     for (i = 0; i < mapsize; i++) {
    256:       PetscReal r = PetscPowReal((PetscReal)R[i] / 255, gamma);
    257:       PetscReal g = PetscPowReal((PetscReal)G[i] / 255, gamma);
    258:       PetscReal b = PetscPowReal((PetscReal)B[i] / 255, gamma);
    259:       R[i]        = (unsigned char)(255 * PetscMin(r, (PetscReal)1.0));
    260:       G[i]        = (unsigned char)(255 * PetscMin(g, (PetscReal)1.0));
    261:       B[i]        = (unsigned char)(255 * PetscMin(b, (PetscReal)1.0));
    262:     }
    263:   }
    264:   PetscFunctionReturn(PETSC_SUCCESS);
    265: }