1 2 #include <../src/sys/classes/draw/utils/axisimpl.h> /*I "petscdraw.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "PetscDrawAxisSetLimits" 6 /*@ 7 PetscDrawAxisSetLimits - Sets the limits (in user coords) of the axis 8 9 Logically Collective on 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 PetscValidHeaderSpecific(axis,PETSC_DRAWAXIS_CLASSID,1); 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)->options,((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 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 PetscValidHeaderSpecific(axis,PETSC_DRAWAXIS_CLASSID,1); 60 *xmin = axis->xlow; 61 *xmax = axis->xhigh; 62 *ymin = axis->ylow; 63 *ymax = axis->yhigh; 64 PetscFunctionReturn(0); 65 } 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "PetscADefLabel" 69 /* 70 val is the label value. sep is the separation to the next (or previous) 71 label; this is useful in determining how many significant figures to 72 keep. 73 */ 74 PetscErrorCode PetscADefLabel(PetscReal val,PetscReal sep,char **p) 75 { 76 static char buf[40]; 77 PetscErrorCode ierr; 78 79 PetscFunctionBegin; 80 /* Find the string */ 81 if (PetscAbsReal(val)/sep < 1.e-4) { 82 buf[0] = '0'; buf[1] = 0; 83 } else { 84 sprintf(buf,"%0.1e",(double)val); 85 ierr = PetscStripZerosPlus(buf);CHKERRQ(ierr); 86 ierr = PetscStripe0(buf);CHKERRQ(ierr); 87 ierr = PetscStripInitialZero(buf);CHKERRQ(ierr); 88 ierr = PetscStripAllZeros(buf);CHKERRQ(ierr); 89 ierr = PetscStripTrailingZeros(buf);CHKERRQ(ierr); 90 } 91 *p =buf; 92 PetscFunctionReturn(0); 93 } 94 95 #undef __FUNCT__ 96 #define __FUNCT__ "PetscADefTicks" 97 /* Finds "nice" locations for the ticks */ 98 PetscErrorCode PetscADefTicks(PetscReal low,PetscReal high,int num,int *ntick,PetscReal *tickloc,int maxtick) 99 { 100 PetscErrorCode ierr; 101 int i,power; 102 PetscReal x = 0.0,base=0.0,eps; 103 104 PetscFunctionBegin; 105 ierr = PetscAGetBase(low,high,num,&base,&power);CHKERRQ(ierr); 106 ierr = PetscAGetNice(low,base,-1,&x);CHKERRQ(ierr); 107 108 /* Values are of the form j * base */ 109 /* Find the starting value */ 110 if (x < low) x += base; 111 112 i = 0; eps = base/10; 113 while (i < maxtick && x <= high+eps) { 114 tickloc[i++] = x; x += base; 115 } 116 *ntick = i; 117 tickloc[i-1] = PetscMin(tickloc[i-1],high); 118 119 if (i < 2 && num < 10) { 120 ierr = PetscADefTicks(low,high,num+1,ntick,tickloc,maxtick);CHKERRQ(ierr); 121 } 122 PetscFunctionReturn(0); 123 } 124 125 #define EPS 1.e-6 126 127 #undef __FUNCT__ 128 #define __FUNCT__ "PetscExp10" 129 PetscErrorCode PetscExp10(PetscReal d,PetscReal *result) 130 { 131 PetscFunctionBegin; 132 *result = PetscPowReal((PetscReal)10.0,d); 133 PetscFunctionReturn(0); 134 } 135 136 #undef __FUNCT__ 137 #define __FUNCT__ "PetscMod" 138 PetscErrorCode PetscMod(PetscReal x,PetscReal y,PetscReal *result) 139 { 140 int i; 141 142 PetscFunctionBegin; 143 if (y == 1) { 144 *result = 0.0; 145 PetscFunctionReturn(0); 146 } 147 i = ((int)x) / ((int)y); 148 x = x - i * y; 149 while (x > y) x -= y; 150 *result = x; 151 PetscFunctionReturn(0); 152 } 153 154 #undef __FUNCT__ 155 #define __FUNCT__ "PetscCopysign" 156 PetscErrorCode PetscCopysign(PetscReal a,PetscReal b,PetscReal *result) 157 { 158 PetscFunctionBegin; 159 if (b >= 0) *result = a; 160 else *result = -a; 161 PetscFunctionReturn(0); 162 } 163 164 #undef __FUNCT__ 165 #define __FUNCT__ "PetscAGetNice" 166 /* 167 Given a value "in" and a "base", return a nice value. 168 based on "sign", extend up (+1) or down (-1) 169 */ 170 PetscErrorCode PetscAGetNice(PetscReal in,PetscReal base,int sign,PetscReal *result) 171 { 172 PetscReal etmp,s,s2,m; 173 PetscErrorCode ierr; 174 175 PetscFunctionBegin; 176 ierr = PetscCopysign (0.5,(double)sign,&s);CHKERRQ(ierr); 177 etmp = in / base + 0.5 + s; 178 ierr = PetscCopysign (0.5,etmp,&s);CHKERRQ(ierr); 179 ierr = PetscCopysign (EPS * etmp,(double)sign,&s2);CHKERRQ(ierr); 180 etmp = etmp - 0.5 + s - s2; 181 ierr = PetscMod(etmp,1.0,&m);CHKERRQ(ierr); 182 etmp = base * (etmp - m); 183 *result = etmp; 184 PetscFunctionReturn(0); 185 } 186 187 #undef __FUNCT__ 188 #define __FUNCT__ "PetscAGetBase" 189 PetscErrorCode PetscAGetBase(PetscReal vmin,PetscReal vmax,int num,PetscReal *Base,int *power) 190 { 191 PetscReal base,ftemp,e10; 192 static PetscReal base_try[5] = {10.0,5.0,2.0,1.0,0.5}; 193 PetscErrorCode ierr; 194 int i; 195 196 PetscFunctionBegin; 197 /* labels of the form n * BASE */ 198 /* get an approximate value for BASE */ 199 base = (vmax - vmin) / (double)(num + 1); 200 201 /* make it of form m x 10^power, m in [1.0, 10) */ 202 if (base <= 0.0) { 203 base = PetscAbsReal(vmin); 204 if (base < 1.0) base = 1.0; 205 } 206 ftemp = PetscLog10Real((1.0 + EPS) * base); 207 if (ftemp < 0.0) ftemp -= 1.0; 208 *power = (int)ftemp; 209 ierr = PetscExp10((double)-*power,&e10);CHKERRQ(ierr); 210 base = base * e10; 211 if (base < 1.0) base = 1.0; 212 /* now reduce it to one of 1, 2, or 5 */ 213 for (i=1; i<5; i++) { 214 if (base >= base_try[i]) { 215 ierr = PetscExp10((double)*power,&e10);CHKERRQ(ierr); 216 base = base_try[i-1] * e10; 217 if (i == 1) *power = *power + 1; 218 break; 219 } 220 } 221 *Base = base; 222 PetscFunctionReturn(0); 223 } 224