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:   Gamma = g;
 19:   return 0;
 20: }

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

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

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

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

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

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

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

 82:   for (i=0; i<mapsize; i++) {
 83:     double u = (double)i/(mapsize-1);
 84:     double m, r=0, g=0, b=0; int k = 0;
 85:     while (k < 4 && u > knots[k+1]) k++;
 86:     m = (u-knots[k])/(knots[k+1]-knots[k]);
 87:     switch(k) {
 88:     case 0: r = 0;     g = 0;   b = (m+1)/2; break;
 89:     case 1: r = 0;     g = m;   b = 1;       break;
 90:     case 2: r = m;     g = 1;   b = 1-m;     break;
 91:     case 3: r = 1;     g = 1-m; b = 0;       break;
 92:     case 4: r = 1-m/2; g = 0;   b = 0;       break;
 93:     }
 94:     R[i] = (unsigned char)(255*PetscMin(r,1.0));
 95:     G[i] = (unsigned char)(255*PetscMin(g,1.0));
 96:     B[i] = (unsigned char)(255*PetscMin(b,1.0));
 97:   }
 98:   return 0;
 99: }

101: static PetscErrorCode PetscDrawCmap_Hot(int mapsize,unsigned char R[],unsigned char G[],unsigned char B[])
102: {
103:   int          i;
104:   const double knots[] =  {0, 3/8., 3/4., 1};

106:   for (i=0; i<mapsize; i++) {
107:     double u = (double)i/(mapsize-1);
108:     double m, r=0, g=0, b=0; int k = 0;
109:     while (k < 2 && u > knots[k+1]) k++;
110:     m = (u-knots[k])/(knots[k+1]-knots[k]);
111:     switch(k) {
112:     case 0: r = m; g = 0; b = 0; break;
113:     case 1: r = 1; g = m; b = 0; break;
114:     case 2: r = 1; g = 1; b = m; break;
115:     }
116:     R[i] = (unsigned char)(255*PetscMin(r,1.0));
117:     G[i] = (unsigned char)(255*PetscMin(g,1.0));
118:     B[i] = (unsigned char)(255*PetscMin(b,1.0));
119:   }
120:   return 0;
121: }

123: static PetscErrorCode PetscDrawCmap_Bone(int mapsize,unsigned char R[],unsigned char G[],unsigned char B[])
124: {
125:   int i;
126:   (void)PetscDrawCmap_Hot(mapsize,R,G,B);
127:   for (i=0; i<mapsize; i++) {
128:     double u = (double)i/(mapsize-1);
129:     double r = (7*u + B[i]/255.0)/8;
130:     double g = (7*u + G[i]/255.0)/8;
131:     double b = (7*u + R[i]/255.0)/8;
132:     R[i] = (unsigned char)(255*PetscMin(r,1.0));
133:     G[i] = (unsigned char)(255*PetscMin(g,1.0));
134:     B[i] = (unsigned char)(255*PetscMin(b,1.0));
135:   }
136:   return 0;
137: }

139: #include "../src/sys/classes/draw/utils/cmap/coolwarm.h"
140: #include "../src/sys/classes/draw/utils/cmap/parula.h"
141: #include "../src/sys/classes/draw/utils/cmap/viridis.h"
142: #include "../src/sys/classes/draw/utils/cmap/plasma.h"
143: #include "../src/sys/classes/draw/utils/cmap/inferno.h"
144: #include "../src/sys/classes/draw/utils/cmap/magma.h"

146: static struct {
147:   const char           *name;
148:   const unsigned char (*data)[3];
149:   PetscErrorCode      (*cmap)(int,unsigned char[],unsigned char[],unsigned char[]);
150: } PetscDrawCmapTable[] = {
151:   {"hue",      NULL, PetscDrawCmap_Hue },     /* varying hue with constant lightness and saturation */
152:   {"gray",     NULL, PetscDrawCmap_Gray},     /* black to white with shades of gray */
153:   {"bone",     NULL, PetscDrawCmap_Bone},     /* black to white with gray-blue shades */
154:   {"jet",      NULL, PetscDrawCmap_Jet },     /* rainbow-like colormap from NCSA, University of Illinois */
155:   {"hot",      NULL, PetscDrawCmap_Hot },     /* black-body radiation */
156:   {"coolwarm", PetscDrawCmap_coolwarm, NULL}, /* ParaView default (Cool To Warm with Diverging interpolation) */
157:   {"parula",   PetscDrawCmap_parula,   NULL}, /* MATLAB (default since R2014b) */
158:   {"viridis",  PetscDrawCmap_viridis,  NULL}, /* matplotlib 1.5 (default since 2.0) */
159:   {"plasma",   PetscDrawCmap_plasma,   NULL}, /* matplotlib 1.5 */
160:   {"inferno",  PetscDrawCmap_inferno,  NULL}, /* matplotlib 1.5 */
161:   {"magma",    PetscDrawCmap_magma,    NULL}, /* matplotlib 1.5 */
162: };

164: PetscErrorCode  PetscDrawUtilitySetCmap(const char colormap[],int mapsize,unsigned char R[],unsigned char G[],unsigned char B[])
165: {
166:   int             i,j;
167:   const char      *cmap_name_list[sizeof(PetscDrawCmapTable)/sizeof(PetscDrawCmapTable[0])];
168:   PetscInt        id = 0, count = (PetscInt)(sizeof(cmap_name_list)/sizeof(char*));
169:   PetscBool       reverse = PETSC_FALSE, brighten = PETSC_FALSE;
170:   PetscReal       beta = 0;

172:   for (i=0; i<count; i++) cmap_name_list[i] = PetscDrawCmapTable[i].name;
173:   if (colormap && colormap[0]) {
174:     PetscBool match = PETSC_FALSE;
175:     for (id=0; !match && id<count; id++) PetscStrcasecmp(colormap,cmap_name_list[id],&match);
177:   }
178:   PetscOptionsGetEList(NULL,NULL,"-draw_cmap",cmap_name_list,count,&id,NULL);
179:   PetscOptionsGetBool(NULL,NULL,"-draw_cmap_reverse",&reverse,NULL);
180:   PetscOptionsGetReal(NULL,NULL,"-draw_cmap_brighten",&beta,&brighten);

183:   if (PetscDrawCmapTable[id].cmap) {
184:     PetscDrawCmapTable[id].cmap(mapsize,R,G,B);
185:   } else {
186:     const unsigned char (*rgb)[3] = PetscDrawCmapTable[id].data;
188:     for (i=0; i<mapsize; i++) {R[i] = rgb[i][0]; G[i] = rgb[i][1]; B[i] = rgb[i][2];}
189:   }

191:   if (reverse) {
192:     i = 0; j = mapsize-1;
193:     while (i < j) {
194: #define SWAP(a,i,j) do { unsigned char t = a[i]; a[i] = a[j]; a[j] = t; } while (0)
195:       SWAP(R,i,j);
196:       SWAP(G,i,j);
197:       SWAP(B,i,j);
198: #undef SWAP
199:       i++; j--;
200:     }
201:   }

203:   if (brighten) {
204:     PetscReal gamma = (beta > 0.0) ? (1 - beta) : (1 / (1 + beta));
205:     for (i=0; i<mapsize; i++) {
206:       PetscReal r = PetscPowReal((PetscReal)R[i]/255,gamma);
207:       PetscReal g = PetscPowReal((PetscReal)G[i]/255,gamma);
208:       PetscReal b = PetscPowReal((PetscReal)B[i]/255,gamma);
209:       R[i] = (unsigned char)(255*PetscMin(r,(PetscReal)1.0));
210:       G[i] = (unsigned char)(255*PetscMin(g,(PetscReal)1.0));
211:       B[i] = (unsigned char)(255*PetscMin(b,(PetscReal)1.0));
212:     }
213:   }
214:   return 0;
215: }