1 #include <petscsys.h> /*I "petscsys.h" I*/ 2 #include <petscdraw.h> 3 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 */ 9 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; 15 16 PetscErrorCode PetscDrawUtilitySetGamma(PetscReal g) { 17 PetscFunctionBegin; 18 Gamma = g; 19 PetscFunctionReturn(0); 20 } 21 22 static inline double PetscHlsHelper(double m1, double m2, double h) { 23 while (h > 1.0) h -= 1.0; 24 while (h < 0.0) h += 1.0; 25 if (h < 1 / 6.0) return m1 + (m2 - m1) * h * 6; 26 if (h < 1 / 2.0) return m2; 27 if (h < 2 / 3.0) return m1 + (m2 - m1) * (2 / 3.0 - h) * 6; 28 return m1; 29 } 30 31 static inline void PetscHlsToRgb(double h, double l, double s, double *r, double *g, double *b) { 32 if (s > 0.0) { 33 double m2 = l <= 0.5 ? l * (1.0 + s) : l + s - (l * s); 34 double m1 = 2 * l - m2; 35 *r = PetscHlsHelper(m1, m2, h + 1 / 3.); 36 *g = PetscHlsHelper(m1, m2, h); 37 *b = PetscHlsHelper(m1, m2, h - 1 / 3.); 38 } else { 39 /* ignore hue */ 40 *r = *g = *b = l; 41 } 42 } 43 44 static inline void PetscGammaCorrect(double *r, double *g, double *b) { 45 PetscReal igamma = 1 / Gamma; 46 *r = (double)PetscPowReal((PetscReal)*r, igamma); 47 *g = (double)PetscPowReal((PetscReal)*g, igamma); 48 *b = (double)PetscPowReal((PetscReal)*b, igamma); 49 } 50 51 static PetscErrorCode PetscDrawCmap_Hue(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 52 int i; 53 double maxhue = 212.0 / 360, lightness = 0.5, saturation = 1.0; 54 55 PetscFunctionBegin; 56 for (i = 0; i < mapsize; i++) { 57 double hue = maxhue * (double)i / (mapsize - 1), r, g, b; 58 PetscHlsToRgb(hue, lightness, saturation, &r, &g, &b); 59 PetscGammaCorrect(&r, &g, &b); 60 R[i] = (unsigned char)(255 * PetscMin(r, 1.0)); 61 G[i] = (unsigned char)(255 * PetscMin(g, 1.0)); 62 B[i] = (unsigned char)(255 * PetscMin(b, 1.0)); 63 } 64 PetscFunctionReturn(0); 65 } 66 67 static PetscErrorCode PetscDrawCmap_Gray(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 68 int i; 69 PetscFunctionBegin; 70 for (i = 0; i < mapsize; i++) R[i] = G[i] = B[i] = (unsigned char)((255.0 * i) / (mapsize - 1)); 71 PetscFunctionReturn(0); 72 } 73 74 static PetscErrorCode PetscDrawCmap_Jet(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 75 int i; 76 const double knots[] = {0, 1 / 8., 3 / 8., 5 / 8., 7 / 8., 1}; 77 78 PetscFunctionBegin; 79 for (i = 0; i < mapsize; i++) { 80 double u = (double)i / (mapsize - 1); 81 double m, r = 0, g = 0, b = 0; 82 int k = 0; 83 while (k < 4 && u > knots[k + 1]) k++; 84 m = (u - knots[k]) / (knots[k + 1] - knots[k]); 85 switch (k) { 86 case 0: 87 r = 0; 88 g = 0; 89 b = (m + 1) / 2; 90 break; 91 case 1: 92 r = 0; 93 g = m; 94 b = 1; 95 break; 96 case 2: 97 r = m; 98 g = 1; 99 b = 1 - m; 100 break; 101 case 3: 102 r = 1; 103 g = 1 - m; 104 b = 0; 105 break; 106 case 4: 107 r = 1 - m / 2; 108 g = 0; 109 b = 0; 110 break; 111 } 112 R[i] = (unsigned char)(255 * PetscMin(r, 1.0)); 113 G[i] = (unsigned char)(255 * PetscMin(g, 1.0)); 114 B[i] = (unsigned char)(255 * PetscMin(b, 1.0)); 115 } 116 PetscFunctionReturn(0); 117 } 118 119 static PetscErrorCode PetscDrawCmap_Hot(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 120 int i; 121 const double knots[] = {0, 3 / 8., 3 / 4., 1}; 122 123 PetscFunctionBegin; 124 for (i = 0; i < mapsize; i++) { 125 double u = (double)i / (mapsize - 1); 126 double m, r = 0, g = 0, b = 0; 127 int k = 0; 128 while (k < 2 && u > knots[k + 1]) k++; 129 m = (u - knots[k]) / (knots[k + 1] - knots[k]); 130 switch (k) { 131 case 0: 132 r = m; 133 g = 0; 134 b = 0; 135 break; 136 case 1: 137 r = 1; 138 g = m; 139 b = 0; 140 break; 141 case 2: 142 r = 1; 143 g = 1; 144 b = m; 145 break; 146 } 147 R[i] = (unsigned char)(255 * PetscMin(r, 1.0)); 148 G[i] = (unsigned char)(255 * PetscMin(g, 1.0)); 149 B[i] = (unsigned char)(255 * PetscMin(b, 1.0)); 150 } 151 PetscFunctionReturn(0); 152 } 153 154 static PetscErrorCode PetscDrawCmap_Bone(int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 155 int i; 156 PetscFunctionBegin; 157 (void)PetscDrawCmap_Hot(mapsize, R, G, B); 158 for (i = 0; i < mapsize; i++) { 159 double u = (double)i / (mapsize - 1); 160 double r = (7 * u + B[i] / 255.0) / 8; 161 double g = (7 * u + G[i] / 255.0) / 8; 162 double b = (7 * u + R[i] / 255.0) / 8; 163 R[i] = (unsigned char)(255 * PetscMin(r, 1.0)); 164 G[i] = (unsigned char)(255 * PetscMin(g, 1.0)); 165 B[i] = (unsigned char)(255 * PetscMin(b, 1.0)); 166 } 167 PetscFunctionReturn(0); 168 } 169 170 #include "../src/sys/classes/draw/utils/cmap/coolwarm.h" 171 #include "../src/sys/classes/draw/utils/cmap/parula.h" 172 #include "../src/sys/classes/draw/utils/cmap/viridis.h" 173 #include "../src/sys/classes/draw/utils/cmap/plasma.h" 174 #include "../src/sys/classes/draw/utils/cmap/inferno.h" 175 #include "../src/sys/classes/draw/utils/cmap/magma.h" 176 177 static struct { 178 const char *name; 179 const unsigned char (*data)[3]; 180 PetscErrorCode (*cmap)(int, unsigned char[], unsigned char[], unsigned char[]); 181 } PetscDrawCmapTable[] = { 182 {"hue", NULL, PetscDrawCmap_Hue }, /* varying hue with constant lightness and saturation */ 183 {"gray", NULL, PetscDrawCmap_Gray}, /* black to white with shades of gray */ 184 {"bone", NULL, PetscDrawCmap_Bone}, /* black to white with gray-blue shades */ 185 {"jet", NULL, PetscDrawCmap_Jet }, /* rainbow-like colormap from NCSA, University of Illinois */ 186 {"hot", NULL, PetscDrawCmap_Hot }, /* black-body radiation */ 187 {"coolwarm", PetscDrawCmap_coolwarm, NULL }, /* ParaView default (Cool To Warm with Diverging interpolation) */ 188 {"parula", PetscDrawCmap_parula, NULL }, /* MATLAB (default since R2014b) */ 189 {"viridis", PetscDrawCmap_viridis, NULL }, /* matplotlib 1.5 (default since 2.0) */ 190 {"plasma", PetscDrawCmap_plasma, NULL }, /* matplotlib 1.5 */ 191 {"inferno", PetscDrawCmap_inferno, NULL }, /* matplotlib 1.5 */ 192 {"magma", PetscDrawCmap_magma, NULL }, /* matplotlib 1.5 */ 193 }; 194 195 PetscErrorCode PetscDrawUtilitySetCmap(const char colormap[], int mapsize, unsigned char R[], unsigned char G[], unsigned char B[]) { 196 int i, j; 197 const char *cmap_name_list[PETSC_STATIC_ARRAY_LENGTH(PetscDrawCmapTable)]; 198 PetscInt id = 0, count = (PetscInt)(sizeof(cmap_name_list) / sizeof(char *)); 199 PetscBool reverse = PETSC_FALSE, brighten = PETSC_FALSE; 200 PetscReal beta = 0; 201 202 PetscFunctionBegin; 203 for (i = 0; i < count; i++) cmap_name_list[i] = PetscDrawCmapTable[i].name; 204 if (colormap && colormap[0]) { 205 PetscBool match = PETSC_FALSE; 206 for (id = 0; !match && id < count; id++) PetscCall(PetscStrcasecmp(colormap, cmap_name_list[id], &match)); 207 PetscCheck(match, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Colormap '%s' not found", colormap); 208 } 209 PetscCall(PetscOptionsGetEList(NULL, NULL, "-draw_cmap", cmap_name_list, count, &id, NULL)); 210 PetscCall(PetscOptionsGetBool(NULL, NULL, "-draw_cmap_reverse", &reverse, NULL)); 211 PetscCall(PetscOptionsGetReal(NULL, NULL, "-draw_cmap_brighten", &beta, &brighten)); 212 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); 213 214 if (PetscDrawCmapTable[id].cmap) { 215 PetscCall(PetscDrawCmapTable[id].cmap(mapsize, R, G, B)); 216 } else { 217 const unsigned char(*rgb)[3] = PetscDrawCmapTable[id].data; 218 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); 219 for (i = 0; i < mapsize; i++) { 220 R[i] = rgb[i][0]; 221 G[i] = rgb[i][1]; 222 B[i] = rgb[i][2]; 223 } 224 } 225 226 if (reverse) { 227 i = 0; 228 j = mapsize - 1; 229 while (i < j) { 230 #define SWAP(a, i, j) \ 231 do { \ 232 unsigned char t = a[i]; \ 233 a[i] = a[j]; \ 234 a[j] = t; \ 235 } while (0) 236 SWAP(R, i, j); 237 SWAP(G, i, j); 238 SWAP(B, i, j); 239 #undef SWAP 240 i++; 241 j--; 242 } 243 } 244 245 if (brighten) { 246 PetscReal gamma = (beta > 0.0) ? (1 - beta) : (1 / (1 + beta)); 247 for (i = 0; i < mapsize; i++) { 248 PetscReal r = PetscPowReal((PetscReal)R[i] / 255, gamma); 249 PetscReal g = PetscPowReal((PetscReal)G[i] / 255, gamma); 250 PetscReal b = PetscPowReal((PetscReal)B[i] / 255, gamma); 251 R[i] = (unsigned char)(255 * PetscMin(r, (PetscReal)1.0)); 252 G[i] = (unsigned char)(255 * PetscMin(g, (PetscReal)1.0)); 253 B[i] = (unsigned char)(255 * PetscMin(b, (PetscReal)1.0)); 254 } 255 } 256 PetscFunctionReturn(0); 257 } 258