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