1 2 #include <../src/sys/classes/draw/utils/axisimpl.h> 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "PetscDrawAxisSetLimits" 6 /*@ 7 PetscDrawAxisSetLimits - Sets the limits (in user coords) of the axis 8 9 Not Collective (ignored on all processors except processor 0 of PetscDrawAxis) 10 11 Input Parameters: 12 + axis - the axis 13 . xmin,xmax - limits in x 14 - ymin,ymax - limits in y 15 16 Options Database: 17 . -drawaxis_hold - hold the initial set of axis limits for future plotting 18 19 Level: advanced 20 21 .seealso: PetscDrawAxisSetHoldLimits() 22 23 @*/ 24 PetscErrorCode PetscDrawAxisSetLimits(PetscDrawAxis axis,PetscReal xmin,PetscReal xmax,PetscReal ymin,PetscReal ymax) 25 { 26 PetscErrorCode ierr; 27 28 PetscFunctionBegin; 29 if (!axis) PetscFunctionReturn(0); 30 if (axis->hold) PetscFunctionReturn(0); 31 axis->xlow = xmin; 32 axis->xhigh= xmax; 33 axis->ylow = ymin; 34 axis->yhigh= ymax; 35 ierr = PetscOptionsHasName(((PetscObject)axis)->prefix,"-drawaxis_hold",&axis->hold);CHKERRQ(ierr); 36 PetscFunctionReturn(0); 37 } 38 39 #undef __FUNCT__ 40 #define __FUNCT__ "PetscDrawAxisGetLimits" 41 /*@ 42 PetscDrawAxisGetLimits - Gets the limits (in user coords) of the axis 43 44 Not Collective (ignored on all processors except processor 0 of PetscDrawAxis) 45 46 Input Parameters: 47 + axis - the axis 48 . xmin,xmax - limits in x 49 - ymin,ymax - limits in y 50 51 Level: advanced 52 53 .seealso: PetscDrawAxisSetLimits() 54 55 @*/ 56 PetscErrorCode PetscDrawAxisGetLimits(PetscDrawAxis axis,PetscReal *xmin,PetscReal *xmax,PetscReal *ymin,PetscReal *ymax) 57 { 58 PetscFunctionBegin; 59 if (!axis) PetscFunctionReturn(0); 60 if (axis->hold) PetscFunctionReturn(0); 61 *xmin = axis->xlow; 62 *xmax = axis->xhigh; 63 *ymin = axis->ylow; 64 *ymax = axis->yhigh; 65 PetscFunctionReturn(0); 66 } 67 68 #undef __FUNCT__ 69 #define __FUNCT__ "PetscADefLabel" 70 /* 71 val is the label value. sep is the separation to the next (or previous) 72 label; this is useful in determining how many significant figures to 73 keep. 74 */ 75 PetscErrorCode PetscADefLabel(PetscReal val,PetscReal sep,char **p) 76 { 77 static char buf[40]; 78 PetscErrorCode ierr; 79 80 PetscFunctionBegin; 81 /* Find the string */ 82 if (PetscAbsReal(val)/sep < 1.e-4) { 83 buf[0] = '0'; buf[1] = 0; 84 } else { 85 sprintf(buf,"%0.1e",(double)val); 86 ierr = PetscStripZerosPlus(buf);CHKERRQ(ierr); 87 ierr = PetscStripe0(buf);CHKERRQ(ierr); 88 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 89 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 90 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 91 } 92 *p =buf; 93 PetscFunctionReturn(0); 94 } 95 96 #undef __FUNCT__ 97 #define __FUNCT__ "PetscADefTicks" 98 /* Finds "nice" locations for the ticks */ 99 PetscErrorCode PetscADefTicks(PetscReal low,PetscReal high,int num,int *ntick,PetscReal * tickloc,int maxtick) 100 { 101 PetscErrorCode ierr; 102 int i,power; 103 PetscReal x = 0.0,base=0.0; 104 105 PetscFunctionBegin; 106 ierr = PetscAGetBase(low,high,num,&base,&power);CHKERRQ(ierr); 107 ierr = PetscAGetNice(low,base,-1,&x);CHKERRQ(ierr); 108 109 /* Values are of the form j * base */ 110 /* Find the starting value */ 111 if (x < low) x += base; 112 113 i = 0; 114 while (i < maxtick && x <= high) { 115 tickloc[i++] = x; 116 x += base; 117 } 118 *ntick = i; 119 120 if (i < 2 && num < 10) { 121 ierr = PetscADefTicks(low,high,num+1,ntick,tickloc,maxtick);CHKERRQ(ierr); 122 } 123 PetscFunctionReturn(0); 124 } 125 126 #define EPS 1.e-6 127 128 #undef __FUNCT__ 129 #define __FUNCT__ "PetscExp10" 130 PetscErrorCode PetscExp10(PetscReal d,PetscReal *result) 131 { 132 PetscFunctionBegin; 133 *result = pow((PetscReal)10.0,d); 134 PetscFunctionReturn(0); 135 } 136 137 #undef __FUNCT__ 138 #define __FUNCT__ "PetscMod" 139 PetscErrorCode PetscMod(PetscReal x,PetscReal y,PetscReal *result) 140 { 141 int i; 142 143 PetscFunctionBegin; 144 if (y == 1) { 145 *result = 0.0; 146 PetscFunctionReturn(0); 147 } 148 i = ((int)x) / ((int)y); 149 x = x - i * y; 150 while (x > y) x -= y; 151 *result = x; 152 PetscFunctionReturn(0); 153 } 154 155 #undef __FUNCT__ 156 #define __FUNCT__ "PetscCopysign" 157 PetscErrorCode PetscCopysign(PetscReal a,PetscReal b,PetscReal *result) 158 { 159 PetscFunctionBegin; 160 if (b >= 0) *result = a; 161 else *result = -a; 162 PetscFunctionReturn(0); 163 } 164 165 #undef __FUNCT__ 166 #define __FUNCT__ "PetscAGetNice" 167 /* 168 Given a value "in" and a "base", return a nice value. 169 based on "sign", extend up (+1) or down (-1) 170 */ 171 PetscErrorCode PetscAGetNice(PetscReal in,PetscReal base,int sign,PetscReal *result) 172 { 173 PetscReal etmp,s,s2,m; 174 PetscErrorCode ierr; 175 176 PetscFunctionBegin; 177 ierr = PetscCopysign (0.5,(double)sign,&s);CHKERRQ(ierr); 178 etmp = in / base + 0.5 + s; 179 ierr = PetscCopysign (0.5,etmp,&s);CHKERRQ(ierr); 180 ierr = PetscCopysign (EPS * etmp,(double)sign,&s2);CHKERRQ(ierr); 181 etmp = etmp - 0.5 + s - s2; 182 ierr = PetscMod(etmp,1.0,&m);CHKERRQ(ierr); 183 etmp = base * (etmp - m); 184 *result = etmp; 185 PetscFunctionReturn(0); 186 } 187 188 #undef __FUNCT__ 189 #define __FUNCT__ "PetscAGetBase" 190 PetscErrorCode PetscAGetBase(PetscReal vmin,PetscReal vmax,int num,PetscReal*Base,int*power) 191 { 192 PetscReal base,ftemp,e10; 193 static PetscReal base_try[5] = {10.0,5.0,2.0,1.0,0.5}; 194 PetscErrorCode ierr; 195 int i; 196 197 PetscFunctionBegin; 198 /* labels of the form n * BASE */ 199 /* get an approximate value for BASE */ 200 base = (vmax - vmin) / (double)(num + 1); 201 202 /* make it of form m x 10^power, m in [1.0, 10) */ 203 if (base <= 0.0) { 204 base = PetscAbsReal(vmin); 205 if (base < 1.0) base = 1.0; 206 } 207 ftemp = log10((1.0 + EPS) * base); 208 if (ftemp < 0.0) ftemp -= 1.0; 209 *power = (int)ftemp; 210 ierr = PetscExp10((double)- *power,&e10);CHKERRQ(ierr); 211 base = base * e10; 212 if (base < 1.0) base = 1.0; 213 /* now reduce it to one of 1, 2, or 5 */ 214 for (i=1; i<5; i++) { 215 if (base >= base_try[i]) { 216 ierr = PetscExp10((double)*power,&e10);CHKERRQ(ierr); 217 base = base_try[i-1] * e10; 218 if (i == 1) *power = *power + 1; 219 break; 220 } 221 } 222 *Base = base; 223 PetscFunctionReturn(0); 224 } 225 226 227 228 229 230 231