File: common\calcfrac.c
1 /*
2 CALCFRAC.C contains the high level ("engine") code for calculating the
3 fractal images (well, SOMEBODY had to do it!).
4 Original author Tim Wegner, but just about ALL the authors have contributed
5 SOME code to this routine at one time or another, or contributed to one of
6 the many massive restructurings.
7 The following modules work very closely with CALCFRAC.C:
8 FRACTALS.C the fractal-specific code for escape-time fractals.
9 FRACSUBR.C assorted subroutines belonging mainly to calcfrac.
10 CALCMAND.ASM fast Mandelbrot/Julia integer implementation
11 Additional fractal-specific modules are also invoked from CALCFRAC:
12 LORENZ.C engine level and fractal specific code for attractors.
13 JB.C julibrot logic
14 PARSER.C formula fractals
15 and more
16 -------------------------------------------------------------------- */
17
18 #include <string.h>
19 #include <limits.h>
20 #include <time.h>
21 /* see Fractint.c for a description of the "include" hierarchy */
22 #include "port.h"
23 #include "prototyp.h"
24 #include "fractype.h"
25 #include "targa_lc.h"
26
27
28 /* routines in this module */
29 static void perform_worklist(void);
30 static int OneOrTwoPass(void);
31 static int _fastcall StandardCalc(int);
32 static int _fastcall potential(double,long);
33 static void decomposition(void);
34 static int bound_trace_main(void);
35 static void step_col_row(void);
36 static int solidguess(void);
37 static int _fastcall guessrow(int,int,int);
38 static void _fastcall plotblock(int,int,int,int);
39 static void _fastcall setsymmetry(int,int);
40 static int _fastcall xsym_split(int,int);
41 static int _fastcall ysym_split(int,int);
42 static void _fastcall puttruecolor_disk(int,int,int);
43 static int diffusion_engine (void);
44 static int sticky_orbits(void);
45
46 /**CJLT new function prototypes: */
47 static int tesseral(void);
48 static int _fastcall tesschkcol(int,int,int);
49 static int _fastcall tesschkrow(int,int,int);
50 static int _fastcall tesscol(int,int,int);
51 static int _fastcall tessrow(int,int,int);
52
53 /* new drawing method by HB */
54 static int diffusion_scan(void);
55
56 /* lookup tables to avoid too much bit fiddling : */
57 char far dif_la[] = {
58 0, 8, 0, 8,4,12,4,12,0, 8, 0, 8,4,12,4,12, 2,10, 2,10,6,14,6,14,2,10,
59 2,10, 6,14,6,14,0, 8,0, 8, 4,12,4,12,0, 8, 0, 8, 4,12,4,12,2,10,2,10,
60 6,14, 6,14,2,10,2,10,6,14, 6,14,1, 9,1, 9, 5,13, 5,13,1, 9,1, 9,5,13,
61 5,13, 3,11,3,11,7,15,7,15, 3,11,3,11,7,15, 7,15, 1, 9,1, 9,5,13,5,13,
62 1, 9, 1, 9,5,13,5,13,3,11, 3,11,7,15,7,15, 3,11, 3,11,7,15,7,15,0, 8,
63 0, 8, 4,12,4,12,0, 8,0, 8, 4,12,4,12,2,10, 2,10, 6,14,6,14,2,10,2,10,
64 6,14, 6,14,0, 8,0, 8,4,12, 4,12,0, 8,0, 8, 4,12, 4,12,2,10,2,10,6,14,
65 6,14, 2,10,2,10,6,14,6,14, 1, 9,1, 9,5,13, 5,13, 1, 9,1, 9,5,13,5,13,
66 3,11, 3,11,7,15,7,15,3,11, 3,11,7,15,7,15, 1, 9, 1, 9,5,13,5,13,1, 9,
67 1, 9, 5,13,5,13,3,11,3,11, 7,15,7,15,3,11, 3,11, 7,15,7,15
68 };
69
70 char far dif_lb[] = {
71 0, 8, 8, 0, 4,12,12, 4, 4,12,12, 4, 8, 0, 0, 8, 2,10,10, 2, 6,14,14,
72 6, 6,14,14, 6,10, 2, 2,10, 2,10,10, 2, 6,14,14, 6, 6,14,14, 6,10, 2,
73 2,10, 4,12,12, 4, 8, 0, 0, 8, 8, 0, 0, 8,12, 4, 4,12, 1, 9, 9, 1, 5,
74 13,13, 5, 5,13,13, 5, 9, 1, 1, 9, 3,11,11, 3, 7,15,15, 7, 7,15,15, 7,
75 11, 3, 3,11, 3,11,11, 3, 7,15,15, 7, 7,15,15, 7,11, 3, 3,11, 5,13,13,
76 5, 9, 1, 1, 9, 9, 1, 1, 9,13, 5, 5,13, 1, 9, 9, 1, 5,13,13, 5, 5,13,
77 13, 5, 9, 1, 1, 9, 3,11,11, 3, 7,15,15, 7, 7,15,15, 7,11, 3, 3,11, 3,
78 11,11, 3, 7,15,15, 7, 7,15,15, 7,11, 3, 3,11, 5,13,13, 5, 9, 1, 1, 9,
79 9, 1, 1, 9,13, 5, 5,13, 2,10,10, 2, 6,14,14, 6, 6,14,14, 6,10, 2, 2,
80 10, 4,12,12, 4, 8, 0, 0, 8, 8, 0, 0, 8,12, 4, 4,12, 4,12,12, 4, 8, 0,
81 0, 8, 8, 0, 0, 8,12, 4, 4,12, 6,14,14, 6,10, 2, 2,10,10, 2, 2,10,14,
82 6, 6,14
83 };
84
85 /* added for testing autologmap() */
86 static long autologmap(void);
87
88 _LCMPLX linitorbit;
89 long lmagnitud, llimit, llimit2, lclosenuff, l16triglim;
90 _CMPLX init,tmp,old,new,saved;
91 int color;
92 long coloriter, oldcoloriter, realcoloriter;
93 int row, col, passes;
94 int iterations, invert;
95 double f_radius,f_xcenter, f_ycenter; /* for inversion */
96 void (_fastcall *putcolor)(int,int,int) = putcolor_a;
97 void (_fastcall *plot)(int,int,int) = putcolor_a;
98 typedef void (_fastcall *PLOTC)(int,int,int);
99 typedef void (_fastcall *GETC)(int,int,int);
100
101 double magnitude, rqlim, rqlim2, rqlim_save;
102 int no_mag_calc = 0;
103 int use_old_period = 0;
104 int use_old_distest = 0;
105 int old_demm_colors = 0;
106 int (*calctype)(void);
107 int (*calctypetmp)(void);
108 int quick_calc = 0;
109 double closeprox = 0.01;
110
111 double closenuff;
112 int pixelpi; /* value of pi in pixels */
113 unsigned long lm; /* magnitude limit (CALCMAND) */
114
115 /* ORBIT variables */
116 int show_orbit; /* flag to turn on and off */
117 int orbit_ptr; /* pointer into save_orbit array */
118 int far *save_orbit; /* array to save orbit values */
119 int orbit_color=15; /* XOR color */
120
121 int ixstart, ixstop, iystart, iystop; /* start, stop here */
122 int symmetry; /* symmetry flag */
123 int reset_periodicity; /* nonzero if escape time pixel rtn to reset */
124 int kbdcount, max_kbdcount; /* avoids checking keyboard too often */
125
126 U16 resume_info = 0; /* handle to resume info if allocated */
127 int resuming; /* nonzero if resuming after interrupt */
128 int num_worklist; /* resume worklist for standard engine */
129 WORKLIST worklist[MAXCALCWORK];
130 int xxstart,xxstop,xxbegin; /* these are same as worklist, */
131 int yystart,yystop,yybegin; /* declared as separate items */
132 int workpass,worksym; /* for the sake of calcmand */
133
134 VOIDFARPTR typespecific_workarea = NULL;
135
136 static double dem_delta, dem_width; /* distance estimator variables */
137 static double dem_toobig;
138 static int dem_mandel;
139 #define DEM_BAILOUT 535.5 /* (pb: not sure if this is special or arbitrary) */
140
141 /* variables which must be visible for tab_display */
142 int got_status; /* -1 if not, 0 for 1or2pass, 1 for ssg, */
143 /* 2 for btm, 3 for 3d, 4 for tesseral, 5 for diffusion_scan */
144 /* 6 for orbits */
145 int curpass,totpasses;
146 int currow,curcol;
147
148 /* static vars for diffusion scan */
149 unsigned bits=0; /* number of bits in the counter */
150 unsigned long dif_counter; /* the diffusion counter */
151 unsigned long dif_limit; /* the diffusion counter */
152
153 /* static vars for solidguess & its subroutines */
154 char three_pass;
155 static int maxblock,halfblock;
156 static int guessplot; /* paint 1st pass row at a time? */
157 static int right_guess,bottom_guess;
158 #define maxyblk 7 /* maxxblk*maxyblk*2 <= 4096, the size of "prefix" */
159 #define maxxblk 202 /* each maxnblk is oversize by 2 for a "border" */
160 /* maxxblk defn must match fracsubr.c */
161 /* next has a skip bit for each maxblock unit;
162 1st pass sets bit [1]... off only if block's contents guessed;
163 at end of 1st pass [0]... bits are set if any surrounding block not guessed;
164 bits are numbered [..][y/16+1][x+1]&(1<<(y&15)) */
165
166 /* Original array */
167 /* extern unsigned int prefix[2][maxyblk][maxxblk]; */
168
169 typedef int (*TPREFIX)[2][maxyblk][maxxblk];
170 #define tprefix (*((TPREFIX)prefix))
171
172 /* size of next puts a limit of MAXPIXELS pixels across on solid guessing logic */
173 #ifdef XFRACT
176 #endif
177
178 int nxtscreenflag; /* for cellular next screen generation */
179 int attractors; /* number of finite attractors */
180 _CMPLX attr[N_ATTR]; /* finite attractor vals (f.p) */
181 _LCMPLX lattr[N_ATTR]; /* finite attractor vals (int) */
182 int attrperiod[N_ATTR]; /* period of the finite attractor */
183
184 /***** vars for new btm *****/
185 enum direction {North,East,South,West};
186 enum direction going_to;
187 int trail_row, trail_col;
188
189 #ifndef sqr
191 #endif
192
193 #ifndef lsqr
195 #endif
196
197 /* -------------------------------------------------------------------- */
198 /* These variables are external for speed's sake only */
199 /* -------------------------------------------------------------------- */
200
201 int periodicitycheck;
202
203 /* For periodicity testing, only in StandardFractal() */
204 int nextsavedincr;
205 long firstsavedand;
206
207 static BYTE *savedots = NULL;
208 static BYTE *fillbuff;
209 static int savedotslen;
210 static int showdotcolor;
211 int atan_colors = 180;
212
213 static int showdot_width = 0;
214 #define SAVE 1
215 #define RESTORE 2
216
217 #define JUST_A_POINT 0
218 #define LOWER_RIGHT 1
219 #define UPPER_RIGHT 2
220 #define LOWER_LEFT 3
221 #define UPPER_LEFT 4
222
223 /* FMODTEST routine. */
224 /* Makes the test condition for the FMOD coloring type
225 that of the current bailout method. 'or' and 'and'
226 methods are not used - in these cases a normal
227 modulus test is used */
228
229 double fmodtest(void)
230 {
231 double result;
232 if (inside==FMODI && save_release <= 2000) /* for backwards compatibility */
233 {
234 if (magnitude == 0.0 || no_mag_calc == 0 || integerfractal)
235 result=sqr(new.x)+sqr(new.y);
236 else
237 result=magnitude; /* don't recalculate */
238 return (result);
239 }
240
241 switch(bailoutest)
242 {
243 case (Mod):
244 {
245 if (magnitude == 0.0 || no_mag_calc == 0 || integerfractal)
246 result=sqr(new.x)+sqr(new.y);
247 else
248 result=magnitude; /* don't recalculate */
249 }break;
250 case (Real):
251 {
252 result=sqr(new.x);
253 }break;
254 case (Imag):
255 {
256 result=sqr(new.y);
257 }break;
258 case (Or):
259 {
260 double tmpx, tmpy;
261 if ((tmpx = sqr(new.x)) > (tmpy = sqr(new.y)))
262 result=tmpx;
263 else
264 result=tmpy;
265 }break;
266 case (Manh):
267 {
268 result=sqr(fabs(new.x)+fabs(new.y));
269 }break;
270 case (Manr):
271 {
272 result=sqr(new.x+new.y);
273 }break;
274 default:
275 {
276 result=sqr(new.x)+sqr(new.y);
277 }break;
278 }
279 return (result);
280 }
281
282 /*
283 The sym_fill_line() routine was pulled out of the boundary tracing
284 code for re-use with showdot. It's purpose is to fill a line with a
285 solid color. This assumes that BYTE *str is already filled
286 with the color. The routine does write the line using symmetry
287 in all cases, however the symmetry logic assumes that the line
288 is one color; it is not general enough to handle a row of
289 pixels of different colors.
290 */
291 static void sym_fill_line(int row, int left, int right, BYTE *str)
292 {
293 int i,j,k, length;
294 length = right-left+1;
295 put_line(row,left,right,str);
296 /* here's where all the symmetry goes */
297 if (plot == putcolor)
298 kbdcount -= length >> 4; /* seems like a reasonable value */
299 else if (plot == symplot2) /* X-axis symmetry */
300 {
301 i = yystop-(row-yystart);
302 if (i > iystop && i < ydots)
303 {
304 put_line(i,left,right,str);
305 kbdcount -= length >> 3;
306 }
307 }
308 else if (plot == symplot2Y) /* Y-axis symmetry */
309 {
310 put_line(row,xxstop-(right-xxstart),xxstop-(left-xxstart),str);
311 kbdcount -= length >> 3;
312 }
313 else if (plot == symplot2J) /* Origin symmetry */
314 {
315 i = yystop-(row-yystart);
316 j = min(xxstop-(right-xxstart),xdots-1);
317 k = min(xxstop-(left -xxstart),xdots-1);
318 if (i > iystop && i < ydots && j <= k)
319 put_line(i,j,k,str);
320 kbdcount -= length >> 3;
321 }
322 else if (plot == symplot4) /* X-axis and Y-axis symmetry */
323 {
324 i = yystop-(row-yystart);
325 j = min(xxstop-(right-xxstart),xdots-1);
326 k = min(xxstop-(left -xxstart),xdots-1);
327 if (i > iystop && i < ydots)
328 {
329 put_line(i,left,right,str);
330 if(j <= k)
331 put_line(i,j,k,str);
332 }
333 if(j <= k)
334 put_line(row,j,k,str);
335 kbdcount -= length >> 2;
336 }
337 else /* cheap and easy way out */
338 {
339 for (i = left; i <= right; i++) /* DG */
340 (*plot)(i,row,str[i-left]);
341 kbdcount -= length >> 1;
342 }
343 }
344
345 /*
346 The sym_put_line() routine is the symmetry-aware version of put_line().
347 It only works efficiently in the no symmetry or XAXIS symmetry case,
348 otherwise it just writes the pixels one-by-one.
349 */
350 static void sym_put_line(int row, int left, int right, BYTE *str)
351 {
352 int length,i;
353 length = right-left+1;
354 put_line(row,left,right,str);
355 if (plot == putcolor)
356 kbdcount -= length >> 4; /* seems like a reasonable value */
357 else if (plot == symplot2) /* X-axis symmetry */
358 {
359 i = yystop-(row-yystart);
360 if (i > iystop && i < ydots)
361 put_line(i,left,right,str);
362 kbdcount -= length >> 3;
363 }
364 else
365 {
366 for (i = left; i <= right; i++) /* DG */
367 (*plot)(i,row,str[i-left]);
368 kbdcount -= length >> 1;
369 }
370 }
371
372 void showdotsaverestore(int startx, int stopx, int starty, int stopy, int direction, int action)
373 {
374 int j,ct;
375 ct = 0;
376 if(direction != JUST_A_POINT)
377 {
378 if(savedots == NULL)
379 {
380 stopmsg(0,"savedots NULL");
381 exit(0);
382 }
383 if(fillbuff == NULL)
384 {
385 stopmsg(0,"fillbuff NULL");
386 exit(0);
387 }
388 }
389 switch(direction)
390 {
391 case LOWER_RIGHT:
392 for(j=starty; j<=stopy; startx++,j++)
393 {
394 if(action==SAVE)
395 {
396 get_line(j,startx,stopx,savedots+ct);
397 sym_fill_line(j,startx,stopx,fillbuff);
398 }
399 else
400 sym_put_line(j,startx,stopx,savedots+ct);
401 ct += stopx-startx+1;
402 }
403 break;
404 case UPPER_RIGHT:
405 for(j=starty; j>=stopy; startx++,j--)
406 {
407 if(action==SAVE)
408 {
409 get_line(j,startx,stopx,savedots+ct);
410 sym_fill_line(j,startx,stopx,fillbuff);
411 }
412 else
413 sym_put_line(j,startx,stopx,savedots+ct);
414 ct += stopx-startx+1;
415 }
416 break;
417 case LOWER_LEFT:
418 for(j=starty; j<=stopy; stopx--,j++)
419 {
420 if(action==SAVE)
421 {
422 get_line(j,startx,stopx,savedots+ct);
423 sym_fill_line(j,startx,stopx,fillbuff);
424 }
425 else
426 sym_put_line(j,startx,stopx,savedots+ct);
427 ct += stopx-startx+1;
428 }
429 break;
430 case UPPER_LEFT:
431 for(j=starty; j>=stopy; stopx--,j--)
432 {
433 if(action==SAVE)
434 {
435 get_line(j,startx,stopx,savedots+ct);
436 sym_fill_line(j,startx,stopx,fillbuff);
437 }
438 else
439 sym_put_line(j,startx,stopx,savedots+ct);
440 ct += stopx-startx+1;
441 }
442 break;
443 }
444 if(action == SAVE)
445 (*plot) (col,row, showdotcolor);
446 }
447
448 int calctypeshowdot(void)
449 {
450 int out, startx, starty, stopx, stopy, direction, width;
451 direction = JUST_A_POINT;
452 startx = stopx = col;
453 starty = stopy = row;
454 width = showdot_width+1;
455 if(width > 0) {
456 if(col+width <= ixstop && row+width <= iystop)
457 {
458 /* preferred showdot shape */
459 direction = UPPER_LEFT;
460 startx = col;
461 stopx = col+width;
462 starty = row+width;
463 stopy = row+1;
464 }
465 else if(col-width >= ixstart && row+width <= iystop)
466 {
467 /* second choice */
468 direction = UPPER_RIGHT;
469 startx = col-width;
470 stopx = col;
471 starty = row+width;
472 stopy = row+1;
473 }
474 else if(col-width >= ixstart && row-width >= iystart)
475 {
476 direction = LOWER_RIGHT;
477 startx = col-width;
478 stopx = col;
479 starty = row-width;
480 stopy = row-1;
481 }
482 else if(col+width <= ixstop && row-width >= iystart)
483 {
484 direction = LOWER_LEFT;
485 startx = col;
486 stopx = col+width;
487 starty = row-width;
488 stopy = row-1;
489 }
490 }
491 showdotsaverestore(startx, stopx, starty, stopy, direction, SAVE);
492 if(orbit_delay > 0) sleepms(orbit_delay);
493 out = (*calctypetmp)();
494 showdotsaverestore(startx, stopx, starty, stopy, direction, RESTORE);
495 return(out);
496 }
497
498 /* use top of extraseg for LogTable if room */
499 int logtable_in_extra_ok(void)
500 {
501 if(((2L*(long)(xdots+ydots)*sizeof(double)+MaxLTSize+1) < (1L<<16))
502 && (bf_math==0)
503 && (fractype != FORMULA)
504 && (fractype != FFORMULA))
505 return(1);
506 else
507 return(0);
508 }
509
510 /******* calcfract - the top level routine for generating an image *******/
511
512 #if (_MSC_VER >= 700)
513 #pragma code_seg ("calcfra1_text") /* place following in an overlay */
514 #endif
515
516 int calcfract(void)
517 {
518 matherr_ct = 0;
519 attractors = 0; /* default to no known finite attractors */
520 display3d = 0;
521 basin = 0;
522 /* added yet another level of indirection to putcolor!!! TW */
523 putcolor = putcolor_a;
524 if (istruecolor && truemode)
525 /* Have to force passes=1 */
526 usr_stdcalcmode = stdcalcmode = '1';
527 if(truecolor)
528 {
529 check_writefile(light_name, ".tga");
530 if(startdisk1(light_name,NULL,0)==0)
531 {
532 /* Have to force passes=1 */
533 usr_stdcalcmode = stdcalcmode = '1';
534 putcolor = puttruecolor_disk;
535 }
536 else
537 truecolor = 0;
538 }
539 if(!use_grid)
540 {
541 if (usr_stdcalcmode != 'o')
542 usr_stdcalcmode = stdcalcmode = '1';
543 }
544
545 init_misc(); /* set up some variables in parser.c */
546 reset_clock();
547
548 /* following delta values useful only for types with rotation disabled */
549 /* currently used only by bifurcation */
550 if (integerfractal)
551 distest = 0;
552 parm.x = param[0];
553 parm.y = param[1];
554 parm2.x = param[2];
555 parm2.y = param[3];
556
557 if (LogFlag && colors < 16) {
558 static FCODE msg[]={"Need at least 16 colors to use logmap"};
559 stopmsg(0,msg);
560 LogFlag = 0;
561 }
562
563 if (use_old_period == 1) {
564 nextsavedincr = 1;
565 firstsavedand = 1;
566 }
567 else {
568 nextsavedincr = (int)log10(maxit); /* works better than log() */
569 if(nextsavedincr < 4) nextsavedincr = 4; /* maintains image with low iterations */
570 firstsavedand = (long)((nextsavedincr*2) + 1);
571 }
572
573 LogTable = NULL;
574 MaxLTSize = maxit;
575 Log_Calc = 0;
576 /* below, INT_MAX=32767 only when an integer is two bytes. Which is not true for Xfractint. */
577 /* Since 32767 is what was meant, replaced the instances of INT_MAX with 32767. */
578 if(LogFlag && (((maxit > 32767) && (save_release > 1920))
579 || Log_Fly_Calc == 1)) {
580 Log_Calc = 1; /* calculate on the fly */
581 SetupLogTable();
582 }
583 else if(LogFlag && (((maxit > 32767) && (save_release <= 1920))
584 || Log_Fly_Calc == 2)) {
585 MaxLTSize = 32767;
586 Log_Calc = 0; /* use logtable */
587 }
588 else if(rangeslen && (maxit >= 32767)) {
589 MaxLTSize = 32766;
590 }
591
592 if ((LogFlag || rangeslen) && !Log_Calc)
593 {
594 if(logtable_in_extra_ok())
595 LogTable = (BYTE far *)(dx0 + 2*(xdots+ydots));
596 else
597 LogTable = (BYTE far *)farmemalloc((long)MaxLTSize + 1);
598
599 if(LogTable == NULL)
600 {
601 if (rangeslen || Log_Fly_Calc == 2) {
602 static FCODE msg[]={"Insufficient memory for logmap/ranges with this maxiter"};
603 stopmsg(0,msg);
604 }
605 else {
606 static FCODE msg[]={"Insufficient memory for logTable, using on-the-fly routine"};
607 stopmsg(0,msg);
608 Log_Fly_Calc = 1;
609 Log_Calc = 1; /* calculate on the fly */
610 SetupLogTable();
611 }
612 }
613 else if (rangeslen) { /* Can't do ranges if MaxLTSize > 32767 */
614 int i,k,l,m,numval,flip,altern;
615 i = k = l = 0;
616 LogFlag = 0; /* ranges overrides logmap */
617 while (i < rangeslen) {
618 m = flip = 0;
619 altern = 32767;
620 if ((numval = ranges[i++]) < 0) {
621 altern = ranges[i++]; /* sub-range iterations */
622 numval = ranges[i++];
623 }
624 if (numval > (int)MaxLTSize || i >= rangeslen)
625 numval = (int)MaxLTSize;
626 while (l <= numval) {
627 LogTable[l++] = (BYTE)(k + flip);
628 if (++m >= altern) {
629 flip ^= 1; /* Alternate colors */
630 m = 0;
631 }
632 }
633 ++k;
634 if (altern != 32767) ++k;
635 }
636 }
637 else
638 SetupLogTable();
639 }
640 lm = 4L << bitshift; /* CALCMAND magnitude limit */
641
642 if (save_release > 2002)
643 atan_colors = colors;
644 else
645 atan_colors = 180;
646
647 /* ORBIT stuff */
648 save_orbit = (int far *)((double huge *)dx0 + 4*OLDMAXPIXELS);
649 show_orbit = start_showorbit;
650 orbit_ptr = 0;
651 orbit_color = 15;
652 if(colors < 16)
653 orbit_color = 1;
654
655 if(inversion[0] != 0.0)
656 {
657 f_radius = inversion[0];
658 f_xcenter = inversion[1];
659 f_ycenter = inversion[2];
660
661 if (inversion[0] == AUTOINVERT) /* auto calc radius 1/6 screen */
662 {
663 inversion[0] = min(fabs(xxmax - xxmin),
664 fabs(yymax - yymin)) / 6.0;
665 fix_inversion(&inversion[0]);
666 f_radius = inversion[0];
667 }
668
669 if (invert < 2 || inversion[1] == AUTOINVERT) /* xcenter not already set */
670 {
671 inversion[1] = (xxmin + xxmax) / 2.0;
672 fix_inversion(&inversion[1]);
673 f_xcenter = inversion[1];
674 if (fabs(f_xcenter) < fabs(xxmax-xxmin) / 100)
675 inversion[1] = f_xcenter = 0.0;
676 }
677
678 if (invert < 3 || inversion[2] == AUTOINVERT) /* ycenter not already set */
679 {
680 inversion[2] = (yymin + yymax) / 2.0;
681 fix_inversion(&inversion[2]);
682 f_ycenter = inversion[2];
683 if (fabs(f_ycenter) < fabs(yymax-yymin) / 100)
684 inversion[2] = f_ycenter = 0.0;
685 }
686
687 invert = 3; /* so values will not be changed if we come back */
688 }
689
690 closenuff = ddelmin*pow(2.0,-(double)(abs(periodicitycheck)));
691 rqlim_save = rqlim;
692 rqlim2 = sqrt(rqlim);
693 if (integerfractal) /* for integer routines (lambda) */
694 {
695 lparm.x = (long)(parm.x * fudge); /* real portion of Lambda */
696 lparm.y = (long)(parm.y * fudge); /* imaginary portion of Lambda */
697 lparm2.x = (long)(parm2.x * fudge); /* real portion of Lambda2 */
698 lparm2.y = (long)(parm2.y * fudge); /* imaginary portion of Lambda2 */
699 llimit = (long)(rqlim * fudge); /* stop if magnitude exceeds this */
700 if (llimit <= 0) llimit = 0x7fffffffL; /* klooge for integer math */
701 llimit2 = (long)(rqlim2 * fudge); /* stop if magnitude exceeds this */
702 lclosenuff = (long)(closenuff * fudge); /* "close enough" value */
703 l16triglim = 8L<<16; /* domain limit of fast trig functions */
704 linitorbit.x = (long)(initorbit.x * fudge);
705 linitorbit.y = (long)(initorbit.y * fudge);
706 }
707 resuming = (calc_status == 2);
708 if (!resuming) /* free resume_info memory if any is hanging around */
709 {
710 end_resume();
711 if (resave_flag) {
712 updatesavename(savename); /* do the pending increment */
713 resave_flag = started_resaves = 0;
714 }
715 calctime = 0;
716 }
717
718 if (curfractalspecific->calctype != StandardFractal
719 && curfractalspecific->calctype != calcmand
720 && curfractalspecific->calctype != calcmandfp
721 && curfractalspecific->calctype != lyapunov
722 && curfractalspecific->calctype != calcfroth)
723 {
724 calctype = curfractalspecific->calctype; /* per_image can override */
725 symmetry = curfractalspecific->symmetry; /* calctype & symmetry */
726 plot = putcolor; /* defaults when setsymmetry not called or does nothing */
727 iystart = ixstart = yystart = xxstart = yybegin = xxbegin = 0;
728 iystop = yystop = ydots -1;
729 ixstop = xxstop = xdots -1;
730 calc_status = 1; /* mark as in-progress */
731 distest = 0; /* only standard escape time engine supports distest */
732 /* per_image routine is run here */
733 if (curfractalspecific->per_image())
734 { /* not a stand-alone */
735 /* next two lines in case periodicity changed */
736 closenuff = ddelmin*pow(2.0,-(double)(abs(periodicitycheck)));
737 lclosenuff = (long)(closenuff * fudge); /* "close enough" value */
738 setsymmetry(symmetry,0);
739 timer(0,calctype); /* non-standard fractal engine */
740 }
741 if (check_key())
742 {
743 if (calc_status == 1) /* calctype didn't set this itself, */
744 calc_status = 3; /* so mark it interrupted, non-resumable */
745 }
746 else
747 calc_status = 4; /* no key, so assume it completed */
748 }
749 else /* standard escape-time engine */
750 {
751 if(stdcalcmode == '3') /* convoluted 'g' + '2' hybrid */
752 {
753 int oldcalcmode;
754 oldcalcmode = stdcalcmode;
755 if(!resuming || three_pass)
756 {
757 stdcalcmode = 'g';
758 three_pass = 1;
759 timer(0,(int(*)())perform_worklist);
760 if(calc_status == 4)
761 {
762 if(xdots >= 640) /* '2' is silly after 'g' for low rez */
763 stdcalcmode = '2';
764 else
765 stdcalcmode = '1';
766
767 timer(0,(int(*)())perform_worklist);
768 three_pass = 0;
769 }
770 }
771 else /* resuming '2' pass */
772 {
773 if(xdots >= 640)
774 stdcalcmode = '2';
775 else
776 stdcalcmode = '1';
777 timer(0,(int(*)())perform_worklist);
778 }
779 stdcalcmode = (char)oldcalcmode;
780 }
781 else /* main case, much nicer! */
782 {
783 three_pass = 0;
784 timer(0,(int(*)())perform_worklist);
785 }
786 }
787 calctime += timer_interval;
788
789 if(LogTable && !Log_Calc)
790 {
791 if(!logtable_in_extra_ok())
792 farmemfree(LogTable); /* free if not using extraseg */
793 LogTable = NULL;
794 }
795 if(typespecific_workarea)
796 {
797 free_workarea();
798 }
799
800 if (curfractalspecific->calctype == calcfroth)
801 froth_cleanup();
802 if((soundflag&7)>1) /* close sound write file */
803 close_snd();
804 if(truecolor)
805 enddisk();
806 return((calc_status == 4) ? 0 : -1);
807 }
808
809 /* locate alternate math record */
810 int find_alternate_math(int type, int math)
811 {
812 int i,ret,curtype /* ,curmath=0 */;
813 /* unsigned umath; */
814 ret = -1;
815 if(math==0)
816 return(ret);
817 i= -1;
818 #if 0 /* for now at least, the only alternatemath is bignum and bigflt math */
830 #else
831 while ((curtype=alternatemath[++i].type) != type && curtype != -1)
832 ;
833 if(curtype == type && alternatemath[i].math)
834 ret = i;
835 #endif
836 return(ret);
837 }
838
839
840 /**************** general escape-time engine routines *********************/
841
842 static void perform_worklist()
843 {
844 int (*sv_orbitcalc)(void) = NULL; /* function that calculates one orbit */
845 int (*sv_per_pixel)(void) = NULL; /* once-per-pixel init */
846 int (*sv_per_image)(void) = NULL; /* once-per-image setup */
847 int i, alt;
848
849 if((alt=find_alternate_math(fractype,bf_math)) > -1)
850 {
851 sv_orbitcalc = curfractalspecific->orbitcalc;
852 sv_per_pixel = curfractalspecific->per_pixel;
853 sv_per_image = curfractalspecific->per_image;
854 curfractalspecific->orbitcalc = alternatemath[alt].orbitcalc;
855 curfractalspecific->per_pixel = alternatemath[alt].per_pixel;
856 curfractalspecific->per_image = alternatemath[alt].per_image;
857 }
858 else
859 bf_math = 0;
860
861 if (potflag && pot16bit)
862 {
863 int tmpcalcmode = stdcalcmode;
864
865 stdcalcmode = '1'; /* force 1 pass */
866 if (resuming == 0)
867 if (pot_startdisk() < 0)
868 {
869 pot16bit = 0; /* startdisk failed or cancelled */
870 stdcalcmode = (char)tmpcalcmode; /* maybe we can carry on??? */
871 }
872 }
873 if (stdcalcmode == 'b' && (curfractalspecific->flags & NOTRACE))
874 stdcalcmode = '1';
875 if (stdcalcmode == 'g' && (curfractalspecific->flags & NOGUESS))
876 stdcalcmode = '1';
877 if (stdcalcmode == 'o' && (curfractalspecific->calctype != StandardFractal))
878 stdcalcmode = '1';
879
880 /* default setup a new worklist */
881 num_worklist = 1;
882 worklist[0].xxstart = worklist[0].xxbegin = 0;
883 worklist[0].yystart = worklist[0].yybegin = 0;
884 worklist[0].xxstop = xdots - 1;
885 worklist[0].yystop = ydots - 1;
886 worklist[0].pass = worklist[0].sym = 0;
887 if (resuming) /* restore worklist, if we can't the above will stay in place */
888 {
889 int vsn;
890 vsn = start_resume();
891 get_resume(sizeof(num_worklist),&num_worklist,sizeof(worklist),worklist,0);
892 end_resume();
893 if (vsn < 2)
894 xxbegin = 0;
895 }
896
897 if (distest) /* setup stuff for distance estimator */
898 {
899 double ftemp,ftemp2,delxx,delyy2,delyy,delxx2,dxsize,dysize;
900 double aspect;
901 if(pseudox && pseudoy)
902 {
903 aspect = (double)pseudoy/(double)pseudox;
904 dxsize = pseudox-1;
905 dysize = pseudoy-1;
906 }
907 else
908 {
909 aspect = (double)ydots/(double)xdots;
910 dxsize = xdots-1;
911 dysize = ydots-1;
912 }
913
914 delxx = (xxmax - xx3rd) / dxsize; /* calculate stepsizes */
915 delyy = (yymax - yy3rd) / dysize;
916 delxx2 = (xx3rd - xxmin) / dysize;
917 delyy2 = (yy3rd - yymin) / dxsize;
918
919 if (save_release < 1827) /* in case it's changed with <G> */
920 use_old_distest = 1;
921 else
922 use_old_distest = 0;
923 rqlim = rqlim_save; /* just in case changed to DEM_BAILOUT earlier */
924 if (distest != 1 || colors == 2) /* not doing regular outside colors */
925 if (rqlim < DEM_BAILOUT) /* so go straight for dem bailout */
926 rqlim = DEM_BAILOUT;
927 if (curfractalspecific->tojulia != NOFRACTAL || use_old_distest
928 || fractype == FORMULA || fractype == FFORMULA)
929 dem_mandel = 1; /* must be mandel type, formula, or old PAR/GIF */
930 else
931 dem_mandel = 0;
932 dem_delta = sqr(delxx) + sqr(delyy2);
933 if ((ftemp = sqr(delyy) + sqr(delxx2)) > dem_delta)
934 dem_delta = ftemp;
935 if (distestwidth == 0)
936 distestwidth = 1;
937 ftemp = distestwidth;
938 if (distestwidth > 0)
939 dem_delta *= sqr(ftemp)/10000; /* multiply by thickness desired */
940 else
941 dem_delta *= 1/(sqr(ftemp)*10000); /* multiply by thickness desired */
942 dem_width = ( sqrt( sqr(xxmax-xxmin) + sqr(xx3rd-xxmin) ) * aspect
943 + sqrt( sqr(yymax-yymin) + sqr(yy3rd-yymin) ) ) / distest;
944 ftemp = (rqlim < DEM_BAILOUT) ? DEM_BAILOUT : rqlim;
945 ftemp += 3; /* bailout plus just a bit */
946 ftemp2 = log(ftemp);
947 if(use_old_distest)
948 dem_toobig = sqr(ftemp) * sqr(ftemp2) * 4 / dem_delta;
949 else
950 dem_toobig = fabs(ftemp) * fabs(ftemp2) * 2 / sqrt(dem_delta);
951 }
952
953 while (num_worklist > 0)
954 {
955 /* per_image can override */
956 calctype = curfractalspecific->calctype;
957 symmetry = curfractalspecific->symmetry; /* calctype & symmetry */
958 plot = putcolor; /* defaults when setsymmetry not called or does nothing */
959
960 /* pull top entry off worklist */
961 ixstart = xxstart = worklist[0].xxstart;
962 ixstop = xxstop = worklist[0].xxstop;
963 xxbegin = worklist[0].xxbegin;
964 iystart = yystart = worklist[0].yystart;
965 iystop = yystop = worklist[0].yystop;
966 yybegin = worklist[0].yybegin;
967 workpass = worklist[0].pass;
968 worksym = worklist[0].sym;
969 --num_worklist;
970 for (i=0; i<num_worklist; ++i)
971 worklist[i] = worklist[i+1];
972
973 calc_status = 1; /* mark as in-progress */
974
975 curfractalspecific->per_image();
976 if(showdot >= 0)
977 {
978 find_special_colors();
979 switch(autoshowdot)
980 {
981 case 'd':
982 showdotcolor = color_dark%colors;
983 break;
984 case 'm':
985 showdotcolor = color_medium%colors;
986 break;
987 case 'b':
988 case 'a':
989 showdotcolor = color_bright%colors;
990 break;
991 default:
992 showdotcolor = showdot%colors;
993 break;
994 }
995 if(sizedot <= 0)
996 showdot_width = -1;
997 else
998 {
999 double dshowdot_width;
1000 dshowdot_width = (double)sizedot*xdots/1024.0;
1001 /*
1002 Arbitrary sanity limit, however showdot_width will
1003 overflow if dshowdot width gets near 256.
1004 */
1005 if(dshowdot_width > 150.0)
1006 showdot_width = 150;
1007 else if(dshowdot_width > 0.0)
1008 showdot_width = (int)dshowdot_width;
1009 else
1010 showdot_width = -1;
1011 }
1012 #ifdef SAVEDOTS_USES_MALLOC
1039 #else
1040 while((savedotslen=sqr(showdot_width)+5*showdot_width+4) > 2048)
1041 showdot_width--;
1042 savedots = (BYTE *)decoderline;
1043 savedotslen /= 2;
1044 fillbuff = savedots + savedotslen;
1045 memset(fillbuff,showdotcolor,savedotslen);
1046 #endif
1047 calctypetmp = calctype;
1048 calctype = calctypeshowdot;
1049 }
1050
1051 /* some common initialization for escape-time pixel level routines */
1052 closenuff = ddelmin*pow(2.0,-(double)(abs(periodicitycheck)));
1053 lclosenuff = (long)(closenuff * fudge); /* "close enough" value */
1054 kbdcount=max_kbdcount;
1055
1056 setsymmetry(symmetry,1);
1057
1058 if (!(resuming)&&(labs(LogFlag) ==2 || (LogFlag && Log_Auto_Calc)))
1059 { /* calculate round screen edges to work out best start for logmap */
1060 LogFlag = ( autologmap() * (LogFlag / labs(LogFlag)));
1061 SetupLogTable();
1062 }
1063
1064 /* call the appropriate escape-time engine */
1065 switch (stdcalcmode)
1066 {
1067 case 's':
1068 if(debugflag==3444)
1069 soi_ldbl();
1070 else
1071 soi();
1072 break;
1073 case 't':
1074 tesseral();
1075 break;
1076 case 'b':
1077 bound_trace_main();
1078 break;
1079 case 'g':
1080 solidguess();
1081 break;
1082 case 'd':
1083 diffusion_scan();
1084 break;
1085 case 'o':
1086 sticky_orbits();
1087 break;
1088 default:
1089 OneOrTwoPass();
1090 }
1091 #ifdef SAVEDOTS_USES_MALLOC
1098 #endif
1099 if (check_key()) /* interrupted? */
1100 break;
1101 }
1102
1103 if (num_worklist > 0)
1104 { /* interrupted, resumable */
1105 alloc_resume(sizeof(worklist)+20,2);
1106 put_resume(sizeof(num_worklist),&num_worklist,sizeof(worklist),worklist,0);
1107 }
1108 else
1109 calc_status = 4; /* completed */
1110 if(sv_orbitcalc != NULL)
1111 {
1112 curfractalspecific->orbitcalc = sv_orbitcalc;
1113 curfractalspecific->per_pixel = sv_per_pixel;
1114 curfractalspecific->per_image = sv_per_image;
1115 }
1116 }
1117
1118 #if (_MSC_VER >= 700)
1119 #pragma code_seg () /* back to normal segment */
1120 #endif
1121
1122 static int diffusion_scan(void)
1123 {
1124 double log2;
1125
1126 log2 = (double) log (2.0);
1127
1128 got_status = 5;
1129
1130 /* note: the max size of 2048x2048 gives us a 22 bit counter that will */
1131 /* fit any 32 bit architecture, the maxinum limit for this case would */
1132 /* be 65536x65536 (HB) */
1133
1134 bits = (unsigned) (min ( log (iystop-iystart+1), log(ixstop-ixstart+1) )/log2 );
1135 bits <<= 1; /* double for two axes */
1136 dif_limit = 1l << bits;
1137
1138 if (diffusion_engine() == -1)
1139 {
1140 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,
1141 (int)(dif_counter >> 16), /* high, */
1142 (int)(dif_counter & 0xffff), /* low order words */
1143 worksym);
1144 return(-1);
1145 }
1146
1147 return(0);
1148 }
1149
1150 /* little macro that plots a filled square of color c, size s with
1151 top left cornet at (x,y) with optimization from sym_fill_line */
1152 #define plot_block(x,y,s,c) \
1153 memset(dstack,(c),(s)); \
1154 for(ty=(y); ty<(y)+(s); ty++) \
1155 sym_fill_line(ty, (x), (x)+(s)-1, dstack)
1156
1157 /* macro that does the same as above, but checks the limits in x and y */
1158 #define plot_block_lim(x,y,s,c) \
1159 memset(dstack,(c),(s)); \
1160 for(ty=(y); ty<min((y)+(s),iystop+1); ty++) \
1161 sym_fill_line(ty, (x), min((x)+(s)-1,ixstop), dstack)
1162
1163 /* macro: count_to_int(dif_counter, colo, rowo) */
1164 /* (inlined function:) */
1165 #define count_to_int(C,x,y)\
1166 \
1167 tC = C; \
1168 \
1169 x =dif_la[tC&0xFF]; y =dif_lb[tC&0xFF]; tC>>=8; \
1170 x <<=4; x+=dif_la[tC&0xFF]; y <<=4; y+=dif_lb[tC&0xFF]; tC>>=8; \
1171 x <<=4; x+=dif_la[tC&0xFF]; y <<=4; y+=dif_lb[tC&0xFF]; tC>>=8; \
1172 x >>= dif_offset; y >>= dif_offset
1173 /* end of inlined function */
1174
1175 /* REMOVED: counter byte 3 */
1176 /* (x) <<=4; (x)+=dif_la[tC&0(x)FF]; (y) <<=4; (y)+=dif_lb[tC&0(x)FF]; tC>>=8;
1177 --> eliminated this and made (*) because fractint user coordinates up to
1178 2048(x)2048 what means a counter of 24 bits or 3 bytes */
1179
1180 /* Calculate the point */
1181 #define calculate \
1182 reset_periodicity = 1; \
1183 if ((*calctype)() == -1) \
1184 return(-1); \
1185 reset_periodicity = 0
1186
1187 static int diffusion_engine (void) {
1188
1189 double log2 = (double) log (2.0);
1190
1191 int i,j;
1192
1193 int nx,ny; /* number of tyles to build in x and y dirs */
1194 /* made this to complete the area that is not */
1195 /* a square with sides like 2 ** n */
1196 int rem_x,rem_y; /* what is left on the last tile to draw */
1197
1198 int ty; /* temp for y */
1199
1200 long unsigned tC; /* temp for dif_counter */
1201
1202 int dif_offset; /* offset for adjusting looked-up values */
1203
1204 int sqsz; /* size of the block being filled */
1205
1206 int colo, rowo; /* original col and row */
1207
1208 int s = 1 << (bits/2); /* size of the square */
1209
1210 nx = (int) floor( (ixstop-ixstart+1)/s );
1211 ny = (int) floor( (iystop-iystart+1)/s );
1212
1213 rem_x = (ixstop-ixstart+1) - nx * s;
1214 rem_y = (iystop-iystart+1) - ny * s;
1215
1216 if (yybegin == iystart && workpass == 0) { /* if restarting on pan: */
1217 dif_counter =0l;
1218 } else {
1219 /* yybegin and passes contain data for resuming the type: */
1220 dif_counter = (((long) ((unsigned)yybegin))<<16) | ((unsigned)workpass);
1221 }
1222
1223 dif_offset = 12-(bits/2); /* offset to adjust coordinates */
1224 /* (*) for 4 bytes use 16 for 3 use 12 etc. */
1225
1226 /*************************************/
1227 /* only the points (dithering only) :*/
1228 if ( fillcolor==0 ){
1229 while (dif_counter < (dif_limit>>1)) {
1230 count_to_int(dif_counter, colo, rowo);
1231
1232 i=0;
1233 col = ixstart + colo; /* get the right tiles */
1234 do {
1235 j=0;
1236 row = iystart + rowo ;
1237 do {
1238 calculate;
1239 (*plot)(col,row,color);
1240 j++;
1241 row += s; /* next tile */
1242 } while (j < ny);
1243 /* in the last y tile we may not need to plot the point
1244 */
1245 if (rowo < rem_y) {
1246 calculate;
1247 (*plot)(col,row,color);
1248 }
1249 i++;
1250 col += s;
1251 } while (i < nx);
1252 /* in the last x tiles we may not need to plot the point */ if
1253 (colo < rem_x) {
1254 row = iystart + rowo;
1255 j=0;
1256 do {
1257 calculate;
1258 (*plot)(col,row,color);
1259 j++;
1260 row += s; /* next tile */
1261 } while (j < ny);
1262 if (rowo < rem_y) {
1263 calculate;
1264 (*plot)(col,row,color);
1265 }
1266 }
1267 dif_counter++;
1268 }
1269 } else {
1270 /*********************************/
1271 /* with progressive filling : */
1272 while (dif_counter < (dif_limit>>1))
1273 {
1274 sqsz=1<<( (int)(bits-(int)( log(dif_counter+0.5)/log2 )-1)/2 );
1275
1276 count_to_int(dif_counter, colo, rowo);
1277
1278 i=0;
1279 do {
1280 j=0;
1281 do {
1282 col = ixstart + colo + i * s; /* get the right tiles */
1283 row = iystart + rowo + j * s;
1284
1285 calculate;
1286 plot_block(col,row,sqsz,color);
1287 j++;
1288 } while (j < ny);
1289 /* in the last tile we may not need to plot the point */
1290 if (rowo < rem_y) {
1291 row = iystart + rowo + ny * s;
1292
1293 calculate;
1294 plot_block_lim(col,row,sqsz,color);
1295 }
1296 i++;
1297 } while (i < nx);
1298 /* in the last tile we may not need to plot the point */
1299 if (colo < rem_x) {
1300 col = ixstart + colo + nx * s;
1301 j=0;
1302 do {
1303 row = iystart + rowo + j * s; /* get the right tiles */
1304
1305 calculate;
1306 plot_block_lim(col,row,sqsz,color);
1307 j++;
1308 } while (j < ny);
1309 if (rowo < rem_y) {
1310 row = iystart + rowo + ny * s;
1311
1312 calculate;
1313 plot_block_lim(col,row,sqsz,color);
1314 }
1315 }
1316
1317 dif_counter++;
1318 }
1319 }
1320 /* from half dif_limit on we only plot 1x1 points :-) */
1321 while (dif_counter < dif_limit)
1322 {
1323 count_to_int(dif_counter, colo, rowo);
1324
1325 i=0;
1326 do {
1327 j=0;
1328 do {
1329 col = ixstart + colo + i * s; /* get the right tiles */
1330 row = iystart + rowo + j * s;
1331
1332 calculate;
1333 (*plot)(col,row,color);
1334 j++;
1335 } while (j < ny);
1336 /* in the last tile we may not need to plot the point */
1337 if (rowo < rem_y) {
1338 row = iystart + rowo + ny * s;
1339
1340 calculate;
1341 (*plot)(col,row,color);
1342 }
1343 i++;
1344 } while (i < nx);
1345 /* in the last tile we may nnt need to plot the point */
1346 if (colo < rem_x) {
1347 col = ixstart + colo + nx * s;
1348 j=0;
1349 do {
1350 row = iystart + rowo + j * s; /* get the right tiles */
1351
1352 calculate;
1353 (*plot)(col,row,color);
1354 j++;
1355 } while (j < ny);
1356 if (rowo < rem_y) {
1357 row = iystart + rowo + ny * s;
1358
1359 calculate;
1360 (*plot)(col,row,color);
1361 }
1362 }
1363 dif_counter++;
1364 }
1365 return(0);
1366 }
1367
1368 /* OLD function (less eficient than the lookup code above:
1369 static void count_to_int (long unsigned C, int *r, int *l) {
1370
1371 int i;
1372
1373 *r = *l = 0;
1374
1375 for (i = bits; i>0; i -= 2){
1376 *r <<=1; *r += C % 2; C >>= 1;
1377 *l <<=1; *l += C % 2; C >>= 1;
1378 }
1379 *l = (*l+*r)%(1<<(bits/2)); * a+b mod 2^k *
1380
1381 }
1382 */
1383
1384 char drawmode = 'r';
1385
1386 static int sticky_orbits(void)
1387 {
1388 got_status = 6; /* for <tab> screen */
1389 totpasses = 1;
1390
1391 if (plotorbits2dsetup() == -1) {
1392 stdcalcmode = 'g';
1393 return(-1);
1394 }
1395
1396 switch (drawmode)
1397 {
1398 case 'r':
1399 default:
1400 /* draw a rectangle */
1401 row = yybegin;
1402 col = xxbegin;
1403
1404 while (row <= iystop)
1405 {
1406 currow = row;
1407 while (col <= ixstop)
1408 {
1409 if (plotorbits2dfloat() == -1)
1410 {
1411 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1412 return(-1); /* interrupted */
1413 }
1414 ++col;
1415 }
1416 col = ixstart;
1417 ++row;
1418 }
1419 break;
1420 case 'l':
1421 {
1422 int dX, dY; /* vector components */
1423 int final, /* final row or column number */
1424 G, /* used to test for new row or column */
1425 inc1, /* G increment when row or column doesn't change */
1426 inc2; /* G increment when row or column changes */
1427 char pos_slope;
1428
1429 dX = ixstop - ixstart; /* find vector components */
1430 dY = iystop - iystart;
1431 pos_slope = (char)(dX > 0); /* is slope positive? */
1432 if (dY < 0)
1433 pos_slope = (char)!pos_slope;
1434 if (abs (dX) > abs (dY)) /* shallow line case */
1435 {
1436 if (dX > 0) /* determine start point and last column */
1437 {
1438 col = xxbegin;
1439 row = yybegin;
1440 final = ixstop;
1441 }
1442 else
1443 {
1444 col = ixstop;
1445 row = iystop;
1446 final = xxbegin;
1447 }
1448 inc1 = 2 * abs (dY); /* determine increments and initial G */
1449 G = inc1 - abs (dX);
1450 inc2 = 2 * (abs (dY) - abs (dX));
1451 if (pos_slope)
1452 while (col <= final) /* step through columns checking for new row */
1453 {
1454 if (plotorbits2dfloat() == -1)
1455 {
1456 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1457 return(-1); /* interrupted */
1458 }
1459 col++;
1460 if (G >= 0) /* it's time to change rows */
1461 {
1462 row++; /* positive slope so increment through the rows */
1463 G += inc2;
1464 }
1465 else /* stay at the same row */
1466 G += inc1;
1467 }
1468 else
1469 while (col <= final) /* step through columns checking for new row */
1470 {
1471 if (plotorbits2dfloat() == -1)
1472 {
1473 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1474 return(-1); /* interrupted */
1475 }
1476 col++;
1477 if (G > 0) /* it's time to change rows */
1478 {
1479 row--; /* negative slope so decrement through the rows */
1480 G += inc2;
1481 }
1482 else /* stay at the same row */
1483 G += inc1;
1484 }
1485 } /* if |dX| > |dY| */
1486 else /* steep line case */
1487 {
1488 if (dY > 0) /* determine start point and last row */
1489 {
1490 col = xxbegin;
1491 row = yybegin;
1492 final = iystop;
1493 }
1494 else
1495 {
1496 col = ixstop;
1497 row = iystop;
1498 final = yybegin;
1499 }
1500 inc1 = 2 * abs (dX); /* determine increments and initial G */
1501 G = inc1 - abs (dY);
1502 inc2 = 2 * (abs (dX) - abs (dY));
1503 if (pos_slope)
1504 while (row <= final) /* step through rows checking for new column */
1505 {
1506 if (plotorbits2dfloat() == -1)
1507 {
1508 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1509 return(-1); /* interrupted */
1510 }
1511 row++;
1512 if (G >= 0) /* it's time to change columns */
1513 {
1514 col++; /* positive slope so increment through the columns */
1515 G += inc2;
1516 }
1517 else /* stay at the same column */
1518 G += inc1;
1519 }
1520 else
1521 while (row <= final) /* step through rows checking for new column */
1522 {
1523 if (plotorbits2dfloat() == -1)
1524 {
1525 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1526 return(-1); /* interrupted */
1527 }
1528 row++;
1529 if (G > 0) /* it's time to change columns */
1530 {
1531 col--; /* negative slope so decrement through the columns */
1532 G += inc2;
1533 }
1534 else /* stay at the same column */
1535 G += inc1;
1536 }
1537 }
1538 } /* end case 'l' */
1539 break;
1540 case 'f': /* this code does not yet work??? */
1541 {
1542 double Xctr,Yctr;
1543 LDBL Magnification; /* LDBL not really needed here, but used to match function parameters */
1544 double Xmagfactor,Rotation,Skew;
1545 int angle;
1546 double factor = PI / 180.0;
1547 double theta;
1548 double xfactor = xdots / 2.0;
1549 double yfactor = ydots / 2.0;
1550
1551 angle = xxbegin; /* save angle in x parameter */
1552
1553 cvtcentermag(&Xctr, &Yctr, &Magnification, &Xmagfactor, &Rotation, &Skew);
1554 if (Rotation <= 0)
1555 Rotation += 360;
1556
1557 while (angle < Rotation)
1558 {
1559 theta = (double)angle * factor;
1560 col = (int)(xfactor + (Xctr + Xmagfactor * cos(theta)));
1561 row = (int)(yfactor + (Yctr + Xmagfactor * sin(theta)));
1562 if (plotorbits2dfloat() == -1)
1563 {
1564 add_worklist(angle,0,0,0,0,0,0,worksym);
1565 return(-1); /* interrupted */
1566 }
1567 angle++;
1568 }
1569
1570
1571 } /* end case 'f' */
1572 break;
1573 } /* end switch */
1574
1575 return(0);
1576 }
1577
1578 static int OneOrTwoPass(void)
1579 {
1580 int i;
1581
1582 totpasses = 1;
1583 if (stdcalcmode == '2') totpasses = 2;
1584 if (stdcalcmode == '2' && workpass == 0) /* do 1st pass of two */
1585 {
1586 if (StandardCalc(1) == -1)
1587 {
1588 add_worklist(xxstart,xxstop,col,yystart,yystop,row,0,worksym);
1589 return(-1);
1590 }
1591 if (num_worklist > 0) /* worklist not empty, defer 2nd pass */
1592 {
1593 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,yystart,1,worksym);
1594 return(0);
1595 }
1596 workpass = 1;
1597 xxbegin = xxstart;
1598 yybegin = yystart;
1599 }
1600 /* second or only pass */
1601 if (StandardCalc(2) == -1)
1602 {
1603 i = yystop;
1604 if (iystop != yystop) /* must be due to symmetry */
1605 i -= row - iystart;
1606 add_worklist(xxstart,xxstop,col,row,i,row,workpass,worksym);
1607 return(-1);
1608 }
1609
1610 return(0);
1611 }
1612
1613 static int _fastcall StandardCalc(int passnum)
1614 {
1615 got_status = 0;
1616 curpass = passnum;
1617 row = yybegin;
1618 col = xxbegin;
1619
1620 while (row <= iystop)
1621 {
1622 currow = row;
1623 reset_periodicity = 1;
1624 while (col <= ixstop)
1625 {
1626 /* on 2nd pass of two, skip even pts */
1627 if (quick_calc && !resuming)
1628 if ((color = getcolor(col,row)) != inside) {
1629 ++col;
1630 continue;
1631 }
1632 if (passnum == 1 || stdcalcmode == '1' || (row&1) != 0 || (col&1) != 0)
1633 {
1634 if ((*calctype)() == -1) /* StandardFractal(), calcmand() or calcmandfp() */
1635 return(-1); /* interrupted */
1636 resuming = 0; /* reset so quick_calc works */
1637 reset_periodicity = 0;
1638 if (passnum == 1) /* first pass, copy pixel and bump col */
1639 {
1640 if ((row&1) == 0 && row < iystop)
1641 {
1642 (*plot)(col,row+1,color);
1643 if ((col&1) == 0 && col < ixstop)
1644 (*plot)(col+1,row+1,color);
1645 }
1646 if ((col&1) == 0 && col < ixstop)
1647 (*plot)(++col,row,color);
1648 }
1649 }
1650 ++col;
1651 }
1652 col = ixstart;
1653 if (passnum == 1 && (row&1) == 0)
1654 ++row;
1655 ++row;
1656 }
1657 return(0);
1658 }
1659
1660
1661 int calcmand(void) /* fast per pixel 1/2/b/g, called with row & col set */
1662 {
1663 /* setup values from far array to avoid using es reg in calcmand.asm */
1664 linitx = lxpixel();
1665 linity = lypixel();
1666 if (calcmandasm() >= 0)
1667 {
1668 if ((LogTable || Log_Calc) /* map color, but not if maxit & adjusted for inside,etc */
1669 && (realcoloriter < maxit || (inside < 0 && coloriter == maxit)))
1670 coloriter = logtablecalc(coloriter);
1671 color = abs((int)coloriter);
1672 if (coloriter >= colors) { /* don't use color 0 unless from inside/outside */
1673 if (save_release <= 1950) {
1674 if (colors < 16)
1675 color &= andcolor;
1676 else
1677 color = ((color - 1) % andcolor) + 1; /* skip color zero */
1678 }
1679 else {
1680 if (colors < 16)
1681 color = (int)(coloriter & andcolor);
1682 else
1683 color = (int)(((coloriter - 1) % andcolor) + 1);
1684 }
1685 }
1686 if(debugflag != 470)
1687 if(color <= 0 && stdcalcmode == 'b' ) /* fix BTM bug */
1688 color = 1;
1689 (*plot) (col, row, color);
1690 }
1691 else
1692 color = (int)coloriter;
1693 return (color);
1694 }
1695
1696 long (*calcmandfpasm)(void);
1697
1698 /************************************************************************/
1699 /* added by Wes Loewer - sort of a floating point version of calcmand() */
1700 /* can also handle invert, any rqlim, potflag, zmag, epsilon cross, */
1701 /* and all the current outside options -Wes Loewer 11/03/91 */
1702 /************************************************************************/
1703 int calcmandfp(void)
1704 {
1705 if(invert)
1706 invertz2(&init);
1707 else
1708 {
1709 init.x = dxpixel();
1710 init.y = dypixel();
1711 }
1712 if (calcmandfpasm() >= 0)
1713 {
1714 if (potflag)
1715 coloriter = potential(magnitude, realcoloriter);
1716 if ((LogTable || Log_Calc) /* map color, but not if maxit & adjusted for inside,etc */
1717 && (realcoloriter < maxit || (inside < 0 && coloriter == maxit)))
1718 coloriter = logtablecalc(coloriter);
1719 color = abs((int)coloriter);
1720 if (coloriter >= colors) { /* don't use color 0 unless from inside/outside */
1721 if (save_release <= 1950) {
1722 if (colors < 16)
1723 color &= andcolor;
1724 else
1725 color = ((color - 1) % andcolor) + 1; /* skip color zero */
1726 }
1727 else {
1728 if (colors < 16)
1729 color = (int)(coloriter & andcolor);
1730 else
1731 color = (int)(((coloriter - 1) % andcolor) + 1);
1732 }
1733 }
1734 if(debugflag != 470)
1735 if(color == 0 && stdcalcmode == 'b' ) /* fix BTM bug */
1736 color = 1;
1737 (*plot) (col, row, color);
1738 }
1739 else
1740 color = (int)coloriter;
1741 return (color);
1742 }
1743 #define STARTRAILMAX FLT_MAX /* just a convenient large number */
1744 #define green 2
1745 #define yellow 6
1746 #if 0
1747 #define NUMSAVED 40 /* define this to save periodicity analysis to file */ 1748 #endif
1749 #if 0
1751 #endif
1752 int StandardFractal(void) /* per pixel 1/2/b/g, called with row & col set */
1753 {
1754 #ifdef NUMSAVED
1759 #endif
1760 long savemaxit;
1761 double tantable[16];
1762 int hooper = 0;
1763 long lcloseprox;
1764 double memvalue = 0.0;
1765 double min_orbit = 100000.0; /* orbit value closest to origin */
1766 long min_index = 0; /* iteration of min_orbit */
1767 long cyclelen = -1;
1768 long savedcoloriter = 0;
1769 int caught_a_cycle;
1770 long savedand;
1771 int savedincr; /* for periodicity checking */
1772 _LCMPLX lsaved;
1773 int i, attracted;
1774 _LCMPLX lat;
1775 _CMPLX at;
1776 _CMPLX deriv;
1777 long dem_color = -1;
1778 _CMPLX dem_new;
1779 int check_freq;
1780 double totaldist = 0.0;
1781 _CMPLX lastz;
1782
1783 lcloseprox = (long)(closeprox*fudge);
1784 savemaxit = maxit;
1785 #ifdef NUMSAVED
1791 #endif
1792 if(inside == STARTRAIL)
1793 {
1794 int i;
1795 for(i=0;i<16;i++)
1796 tantable[i] = 0.0;
1797 if (save_release > 1824) maxit = 16;
1798 }
1799 if (periodicitycheck == 0 || inside == ZMAG || inside == STARTRAIL)
1800 oldcoloriter = 2147483647L; /* don't check periodicity at all */
1801 else if (inside == PERIOD) /* for display-periodicity */
1802 oldcoloriter = (maxit/5)*4; /* don't check until nearly done */
1803 else if (reset_periodicity)
1804 oldcoloriter = 255; /* don't check periodicity 1st 250 iterations */
1805
1806 /* Jonathan - how about this idea ? skips first saved value which never works */
1807 #ifdef MINSAVEDAND
1810 #else
1811 if (oldcoloriter < firstsavedand) /* I like it! */
1812 oldcoloriter = firstsavedand;
1813 #endif
1814 /* really fractal specific, but we'll leave it here */
1815 if (!integerfractal)
1816 {
1817 if (useinitorbit == 1)
1818 saved = initorbit;
1819 else {
1820 saved.x = 0;
1821 saved.y = 0;
1822 }
1823 #ifdef NUMSAVED
1825 #endif
1826 if(bf_math)
1827 {
1828 if(decimals > 200)
1829 kbdcount = -1;
1830 if (bf_math == BIGNUM)
1831 {
1832 clear_bn(bnsaved.x);
1833 clear_bn(bnsaved.y);
1834 }
1835 else if (bf_math == BIGFLT)
1836 {
1837 clear_bf(bfsaved.x);
1838 clear_bf(bfsaved.y);
1839 }
1840 }
1841 init.y = dypixel();
1842 if (distest)
1843 {
1844 if (use_old_distest) {
1845 rqlim = rqlim_save;
1846 if (distest != 1 || colors == 2) /* not doing regular outside colors */
1847 if (rqlim < DEM_BAILOUT) /* so go straight for dem bailout */
1848 rqlim = DEM_BAILOUT;
1849 dem_color = -1;
1850 }
1851 deriv.x = 1;
1852 deriv.y = 0;
1853 magnitude = 0;
1854 }
1855 }
1856 else
1857 {
1858 if (useinitorbit == 1)
1859 lsaved = linitorbit;
1860 else {
1861 lsaved.x = 0;
1862 lsaved.y = 0;
1863 }
1864 linit.y = lypixel();
1865 }
1866 orbit_ptr = 0;
1867 coloriter = 0;
1868 if(fractype==JULIAFP || fractype==JULIA)
1869 coloriter = -1;
1870 caught_a_cycle = 0;
1871 if (inside == PERIOD) {
1872 savedand = 16; /* begin checking every 16th cycle */
1873 } else {
1874 /* Jonathan - don't understand such a low savedand -- how about this? */
1875 #ifdef MINSAVEDAND
1877 #else
1878 savedand = firstsavedand; /* begin checking every other cycle */
1879 #endif
1880 }
1881 savedincr = 1; /* start checking the very first time */
1882
1883 if (inside <= BOF60 && inside >= BOF61)
1884 {
1885 magnitude = lmagnitud = 0;
1886 min_orbit = 100000.0;
1887 }
1888 overflow = 0; /* reset integer math overflow flag */
1889
1890 curfractalspecific->per_pixel(); /* initialize the calculations */
1891
1892 attracted = FALSE;
1893
1894 if (outside == TDIS) {
1895 if(integerfractal)
1896 {
1897 old.x = ((double)lold.x) / fudge;
1898 old.y = ((double)lold.y) / fudge;
1899 }
1900 else if (bf_math == BIGNUM)
1901 old = cmplxbntofloat(&bnold);
1902 else if (bf_math==BIGFLT)
1903 old = cmplxbftofloat(&bfold);
1904 lastz.x = old.x;
1905 lastz.y = old.y;
1906 }
1907
1908 if (((soundflag&7) > 2 || showdot >= 0) && orbit_delay > 0)
1909 check_freq = 16;
1910 else
1911 check_freq = 2048;
1912
1913 if(show_orbit)
1914 snd_time_write();
1915 while (++coloriter < maxit)
1916 {
1917 /* calculation of one orbit goes here */
1918 /* input in "old" -- output in "new" */
1919 if (coloriter % check_freq == 0) {
1920 if (check_key())
1921 return (-1);
1922 }
1923
1924 if (distest)
1925 {
1926 double ftemp;
1927 /* Distance estimator for points near Mandelbrot set */
1928 /* Original code by Phil Wilson, hacked around by PB */
1929 /* Algorithms from Peitgen & Saupe, Science of Fractal Images, p.198 */
1930 if (dem_mandel)
1931 ftemp = 2 * (old.x * deriv.x - old.y * deriv.y) + 1;
1932 else
1933 ftemp = 2 * (old.x * deriv.x - old.y * deriv.y);
1934 deriv.y = 2 * (old.y * deriv.x + old.x * deriv.y);
1935 deriv.x = ftemp;
1936 if (use_old_distest) {
1937 if (sqr(deriv.x)+sqr(deriv.y) > dem_toobig)
1938 break;
1939 }
1940 else if (save_release > 1950)
1941 if (max(fabs(deriv.x),fabs(deriv.y)) > dem_toobig)
1942 break;
1943 /* if above exit taken, the later test vs dem_delta will place this
1944 point on the boundary, because mag(old)<bailout just now */
1945
1946 if (curfractalspecific->orbitcalc() || (overflow && save_release > 1826))
1947 {
1948 if (use_old_distest) {
1949 if (dem_color < 0) {
1950 dem_color = coloriter;
1951 dem_new = new;
1952 }
1953 if (rqlim >= DEM_BAILOUT
1954 || magnitude >= (rqlim = DEM_BAILOUT)
1955 || magnitude == 0)
1956 break;
1957 }
1958 else
1959 break;
1960 }
1961 old = new;
1962 }
1963
1964 /* the usual case */
1965 else if ((curfractalspecific->orbitcalc() && inside != STARTRAIL)
1966 || overflow)
1967 break;
1968 if (show_orbit) {
1969 if (!integerfractal)
1970 {
1971 if (bf_math == BIGNUM)
1972 new = cmplxbntofloat(&bnnew);
1973 else if (bf_math==BIGFLT)
1974 new = cmplxbftofloat(&bfnew);
1975 plot_orbit(new.x, new.y, -1);
1976 }
1977 else
1978 iplot_orbit(lnew.x, lnew.y, -1);
1979 }
1980 if( inside < -1)
1981 {
1982 if (bf_math == BIGNUM)
1983 new = cmplxbntofloat(&bnnew);
1984 else if (bf_math == BIGFLT)
1985 new = cmplxbftofloat(&bfnew);
1986 if(inside == STARTRAIL)
1987 {
1988 if(0 < coloriter && coloriter < 16)
1989 {
1990 if (integerfractal)
1991 {
1992 new.x = lnew.x;
1993 new.x /= fudge;
1994 new.y = lnew.y;
1995 new.y /= fudge;
1996 }
1997
1998 if (save_release > 1824) {
1999 if(new.x > STARTRAILMAX)
2000 new.x = STARTRAILMAX;
2001 if(new.x < -STARTRAILMAX)
2002 new.x = -STARTRAILMAX;
2003 if(new.y > STARTRAILMAX)
2004 new.y = STARTRAILMAX;
2005 if(new.y < -STARTRAILMAX)
2006 new.y = -STARTRAILMAX;
2007 tempsqrx = new.x * new.x;
2008 tempsqry = new.y * new.y;
2009 magnitude = tempsqrx + tempsqry;
2010 old = new;
2011 }
2012 {
2013 int tmpcolor;
2014 tmpcolor = (int)(((coloriter - 1) % andcolor) + 1);
2015 tantable[tmpcolor-1] = new.y/(new.x+.000001);
2016 }
2017 }
2018 }
2019 else if(inside == EPSCROSS)
2020 {
2021 hooper = 0;
2022 if(integerfractal)
2023 {
2024 if(labs(lnew.x) < labs(lcloseprox))
2025 {
2026 hooper = (lcloseprox>0? 1 : -1); /* close to y axis */
2027 goto plot_inside;
2028 }
2029 else if(labs(lnew.y) < labs(lcloseprox))
2030 {
2031 hooper = (lcloseprox>0 ? 2: -2); /* close to x axis */
2032 goto plot_inside;
2033 }
2034 }
2035 else
2036 {
2037 if(fabs(new.x) < fabs(closeprox))
2038 {
2039 hooper = (closeprox>0? 1 : -1); /* close to y axis */
2040 goto plot_inside;
2041 }
2042 else if(fabs(new.y) < fabs(closeprox))
2043 {
2044 hooper = (closeprox>0? 2 : -2); /* close to x axis */
2045 goto plot_inside;
2046 }
2047 }
2048 }
2049 else if (inside == FMODI)
2050 {
2051 double mag;
2052 if(integerfractal)
2053 {
2054 new.x = ((double)lnew.x) / fudge;
2055 new.y = ((double)lnew.y) / fudge;
2056 }
2057 mag = fmodtest();
2058 if(mag < closeprox)
2059 memvalue = mag;
2060 }
2061 else if (inside <= BOF60 && inside >= BOF61)
2062 {
2063 if (integerfractal)
2064 {
2065 if (lmagnitud == 0 || no_mag_calc == 0)
2066 lmagnitud = lsqr(lnew.x) + lsqr(lnew.y);
2067 magnitude = lmagnitud;
2068 magnitude = magnitude / fudge;
2069 }
2070 else
2071 if (magnitude == 0.0 || no_mag_calc == 0)
2072 magnitude = sqr(new.x) + sqr(new.y);
2073 if (magnitude < min_orbit)
2074 {
2075 min_orbit = magnitude;
2076 min_index = coloriter + 1;
2077 }
2078 }
2079 }
2080
2081 if (outside == TDIS || outside == FMOD)
2082 {
2083 if (bf_math == BIGNUM)
2084 new = cmplxbntofloat(&bnnew);
2085 else if (bf_math == BIGFLT)
2086 new = cmplxbftofloat(&bfnew);
2087 if (outside == TDIS)
2088 {
2089 if(integerfractal)
2090 {
2091 new.x = ((double)lnew.x) / fudge;
2092 new.y = ((double)lnew.y) / fudge;
2093 }
2094 totaldist += sqrt(sqr(lastz.x-new.x)+sqr(lastz.y-new.y));
2095 lastz.x = new.x;
2096 lastz.y = new.y;
2097 }
2098 else if (outside == FMOD)
2099 {
2100 double mag;
2101 if(integerfractal)
2102 {
2103 new.x = ((double)lnew.x) / fudge;
2104 new.y = ((double)lnew.y) / fudge;
2105 }
2106 mag = fmodtest();
2107 if(mag < closeprox)
2108 memvalue = mag;
2109 }
2110 }
2111
2112 if (attractors > 0) /* finite attractor in the list */
2113 { /* NOTE: Integer code is UNTESTED */
2114 if (integerfractal)
2115 {
2116 for (i = 0; i < attractors; i++)
2117 {
2118 lat.x = lnew.x - lattr[i].x;
2119 lat.x = lsqr(lat.x);
2120 if (lat.x < l_at_rad)
2121 {
2122 lat.y = lnew.y - lattr[i].y;
2123 lat.y = lsqr(lat.y);
2124 if (lat.y < l_at_rad)
2125 {
2126 if ((lat.x + lat.y) < l_at_rad)
2127 {
2128 attracted = TRUE;
2129 if (finattract<0) coloriter = (coloriter%attrperiod[i])+1;
2130 break;
2131 }
2132 }
2133 }
2134 }
2135 }
2136 else
2137 {
2138 for (i = 0; i < attractors; i++)
2139 {
2140 at.x = new.x - attr[i].x;
2141 at.x = sqr(at.x);
2142 if (at.x < f_at_rad)
2143 {
2144 at.y = new.y - attr[i].y;
2145 at.y = sqr(at.y);
2146 if ( at.y < f_at_rad)
2147 {
2148 if ((at.x + at.y) < f_at_rad)
2149 {
2150 attracted = TRUE;
2151 if (finattract<0) coloriter = (coloriter%attrperiod[i])+1;
2152 break;
2153 }
2154 }
2155 }
2156 }
2157 }
2158 if (attracted)
2159 break; /* AHA! Eaten by an attractor */
2160 }
2161
2162 if (coloriter > oldcoloriter) /* check periodicity */
2163 {
2164 if ((coloriter & savedand) == 0) /* time to save a new value */
2165 {
2166 savedcoloriter = coloriter;
2167 if (integerfractal)
2168 lsaved = lnew;/* integer fractals */
2169 else if (bf_math == BIGNUM)
2170 {
2171 copy_bn(bnsaved.x,bnnew.x);
2172 copy_bn(bnsaved.y,bnnew.y);
2173 }
2174 else if (bf_math == BIGFLT)
2175 {
2176 copy_bf(bfsaved.x,bfnew.x);
2177 copy_bf(bfsaved.y,bfnew.y);
2178 }
2179 else
2180 {
2181 saved = new; /* floating pt fractals */
2182 #ifdef NUMSAVED
2188 #endif
2189 }
2190 if (--savedincr == 0) /* time to lengthen the periodicity? */
2191 {
2192 savedand = (savedand << 1) + 1; /* longer periodicity */
2193 savedincr = nextsavedincr;/* restart counter */
2194 }
2195 }
2196 else /* check against an old save */
2197 {
2198 if (integerfractal) /* floating-pt periodicity chk */
2199 {
2200 if (labs(lsaved.x - lnew.x) < lclosenuff)
2201 if (labs(lsaved.y - lnew.y) < lclosenuff)
2202 caught_a_cycle = 1;
2203 }
2204 else if (bf_math == BIGNUM)
2205 {
2206 if (cmp_bn(abs_a_bn(sub_bn(bntmp,bnsaved.x,bnnew.x)), bnclosenuff) < 0)
2207 if (cmp_bn(abs_a_bn(sub_bn(bntmp,bnsaved.y,bnnew.y)), bnclosenuff) < 0)
2208 caught_a_cycle = 1;
2209 }
2210 else if (bf_math == BIGFLT)
2211 {
2212 if (cmp_bf(abs_a_bf(sub_bf(bftmp,bfsaved.x,bfnew.x)), bfclosenuff) < 0)
2213 if (cmp_bf(abs_a_bf(sub_bf(bftmp,bfsaved.y,bfnew.y)), bfclosenuff) < 0)
2214 caught_a_cycle = 1;
2215 }
2216 else
2217 {
2218 if (fabs(saved.x - new.x) < closenuff)
2219 if (fabs(saved.y - new.y) < closenuff)
2220 caught_a_cycle = 1;
2221 #ifdef NUMSAVED
2232 #endif
2233 }
2234 if(caught_a_cycle)
2235 {
2236 #ifdef NUMSAVED
2242 #endif
2243 cyclelen = coloriter-savedcoloriter;
2244 #ifdef NUMSAVED
2245 fprintf(fp,"row %3d col %3d len %6ld iter %6ld savedand %6ld\n",
2246 row,col,cyclelen,coloriter,savedand);
2247 if(zctr > 1 && zctr < NUMSAVED)
2248 {
2249 int i;
2250 for(i=0;i<zctr;i++)
2251 fprintf(fp," caught %2d saved %6ld iter %6ld\n",i,changed[i],caught[i]);
2252 }
2253 fflush(fp); 2254 #endif
2255 coloriter = maxit - 1;
2256 }
2257
2258 }
2259 }
2260 } /* end while (coloriter++ < maxit) */
2261
2262 if (show_orbit)
2263 scrub_orbit();
2264
2265 realcoloriter = coloriter; /* save this before we start adjusting it */
2266 if (coloriter >= maxit)
2267 oldcoloriter = 0; /* check periodicity immediately next time */
2268 else
2269 {
2270 oldcoloriter = coloriter + 10; /* check when past this + 10 next time */
2271 if (coloriter == 0)
2272 coloriter = 1; /* needed to make same as calcmand */
2273 }
2274
2275 if (potflag)
2276 {
2277 if (integerfractal) /* adjust integer fractals */
2278 {
2279 new.x = ((double)lnew.x) / fudge;
2280 new.y = ((double)lnew.y) / fudge;
2281 }
2282 else if (bf_math==BIGNUM)
2283 {
2284 new.x = (double)bntofloat(bnnew.x);
2285 new.y = (double)bntofloat(bnnew.y);
2286 }
2287 else if (bf_math==BIGFLT)
2288 {
2289 new.x = (double)bftofloat(bfnew.x);
2290 new.y = (double)bftofloat(bfnew.y);
2291 }
2292 magnitude = sqr(new.x) + sqr(new.y);
2293 coloriter = potential(magnitude, coloriter);
2294 if (LogTable || Log_Calc)
2295 coloriter = logtablecalc(coloriter);
2296 goto plot_pixel; /* skip any other adjustments */
2297 }
2298
2299 if (coloriter >= maxit) /* an "inside" point */
2300 goto plot_inside; /* distest, decomp, biomorph don't apply */
2301
2302
2303 if (outside < -1) /* these options by Richard Hughes modified by TW */
2304 {
2305 if (integerfractal)
2306 {
2307 new.x = ((double)lnew.x) / fudge;
2308 new.y = ((double)lnew.y) / fudge;
2309 }
2310 else if(bf_math==1)
2311 {
2312 new.x = (double)bntofloat(bnnew.x);
2313 new.y = (double)bntofloat(bnnew.y);
2314 }
2315 /* Add 7 to overcome negative values on the MANDEL */
2316 if (outside == REAL) /* "real" */
2317 coloriter += (long)new.x + 7;
2318 else if (outside == IMAG) /* "imag" */
2319 coloriter += (long)new.y + 7;
2320 else if (outside == MULT && new.y) /* "mult" */
2321 coloriter = (long)((double)coloriter * (new.x/new.y));
2322 else if (outside == SUM) /* "sum" */
2323 coloriter += (long)(new.x + new.y);
2324 else if (outside == ATAN) /* "atan" */
2325 coloriter = (long)fabs(atan2(new.y,new.x)*atan_colors/PI);
2326 else if (outside == FMOD)
2327 coloriter = (long)(memvalue * colors / closeprox);
2328 else if (outside == TDIS) {
2329 coloriter = (long)(totaldist);
2330 }
2331
2332
2333 /* eliminate negative colors & wrap arounds */
2334 if ((coloriter <= 0 || coloriter > maxit) && outside!=FMOD) {
2335 if (save_release < 1961)
2336 coloriter = 0;
2337 else
2338 coloriter = 1;
2339 }
2340 }
2341
2342 if (distest)
2343 {
2344 double dist,temp;
2345 dist = sqr(new.x) + sqr(new.y);
2346 if (dist == 0 || overflow)
2347 dist = 0;
2348 else {
2349 temp = log(dist);
2350 dist = dist * sqr(temp) / ( sqr(deriv.x) + sqr(deriv.y) );
2351 }
2352 if (dist < dem_delta) /* point is on the edge */
2353 {
2354 if (distest > 0)
2355 goto plot_inside; /* show it as an inside point */
2356 coloriter = 0 - distest; /* show boundary as specified color */
2357 goto plot_pixel; /* no further adjustments apply */
2358 }
2359 if (colors == 2)
2360 {
2361 coloriter = !inside; /* the only useful distest 2 color use */
2362 goto plot_pixel; /* no further adjustments apply */
2363 }
2364 if (distest > 1) /* pick color based on distance */
2365 {
2366 if (old_demm_colors) /* this one is needed for old color scheme */
2367 coloriter = (long)sqrt(sqrt(dist) / dem_width + 1);
2368 else if (use_old_distest)
2369 coloriter = (long)sqrt(dist / dem_width + 1);
2370 else
2371 coloriter = (long)(dist / dem_width + 1);
2372 coloriter &= LONG_MAX; /* oops - color can be negative */
2373 goto plot_pixel; /* no further adjustments apply */
2374 }
2375 if (use_old_distest) {
2376 coloriter = dem_color;
2377 new = dem_new;
2378 }
2379 /* use pixel's "regular" color */
2380 }
2381
2382 if (decomp[0] > 0)
2383 decomposition();
2384 else if (biomorph != -1)
2385 {
2386 if (integerfractal)
2387 {
2388 if (labs(lnew.x) < llimit2 || labs(lnew.y) < llimit2)
2389 coloriter = biomorph;
2390 }
2391 else
2392 if (fabs(new.x) < rqlim2 || fabs(new.y) < rqlim2)
2393 coloriter = biomorph;
2394 }
2395
2396 if (outside >= 0 && attracted == FALSE) /* merge escape-time stripes */
2397 coloriter = outside;
2398 else if (LogTable || Log_Calc)
2399 coloriter = logtablecalc(coloriter);
2400 goto plot_pixel;
2401
2402 plot_inside: /* we're "inside" */
2403 if (periodicitycheck < 0 && caught_a_cycle)
2404 coloriter = 7; /* show periodicity */
2405 else if (inside >= 0)
2406 coloriter = inside; /* set to specified color, ignore logpal */
2407 else
2408 {
2409 if(inside == STARTRAIL)
2410 {
2411 int i;
2412 double diff;
2413 coloriter = 0;
2414 for(i=1;i<16;i++)
2415 {
2416 diff = tantable[0] - tantable[i];
2417 if(fabs(diff) < .05)
2418 {
2419 coloriter = i;
2420 break;
2421 }
2422 }
2423 }
2424 else if(inside== PERIOD) {
2425 if (cyclelen>0) {
2426 coloriter = cyclelen;
2427 } else {
2428 coloriter = maxit;
2429 }
2430 }
2431 else if(inside == EPSCROSS)
2432 {
2433 if(hooper==1)
2434 coloriter = green;
2435 else if(hooper==2)
2436 coloriter = yellow;
2437 else if( hooper==0)
2438 coloriter = maxit;
2439 if (show_orbit)
2440 scrub_orbit();
2441 }
2442 else if(inside == FMODI)
2443 {
2444 coloriter = (long)(memvalue * colors / closeprox);
2445 }
2446 else if (inside == ATANI) /* "atan" */
2447 if (integerfractal) {
2448 new.x = ((double)lnew.x) / fudge;
2449 new.y = ((double)lnew.y) / fudge;
2450 coloriter = (long)fabs(atan2(new.y,new.x)*atan_colors/PI);
2451 }
2452 else
2453 coloriter = (long)fabs(atan2(new.y,new.x)*atan_colors/PI);
2454 else if (inside == BOF60)
2455 coloriter = (long)(sqrt(min_orbit) * 75);
2456 else if (inside == BOF61)
2457 coloriter = min_index;
2458 else if (inside == ZMAG)
2459 {
2460 if (integerfractal)
2461 {
2462 /*
2463 new.x = ((double)lnew.x) / fudge;
2464 new.y = ((double)lnew.y) / fudge;
2465 coloriter = (long)((((double)lsqr(lnew.x))/fudge + ((double)lsqr(lnew.y))/fudge) * (maxit>>1) + 1);
2466 */
2467 coloriter = (long)(((double)lmagnitud/fudge) * (maxit>>1) + 1);
2468 }
2469 else
2470 coloriter = (long)((sqr(new.x) + sqr(new.y)) * (maxit>>1) + 1);
2471 }
2472 else /* inside == -1 */
2473 coloriter = maxit;
2474 if (LogTable || Log_Calc)
2475 coloriter = logtablecalc(coloriter);
2476 }
2477
2478 plot_pixel:
2479
2480 color = abs((int)coloriter);
2481 if (coloriter >= colors) { /* don't use color 0 unless from inside/outside */
2482 if (save_release <= 1950) {
2483 if (colors < 16)
2484 color &= andcolor;
2485 else
2486 color = ((color - 1) % andcolor) + 1; /* skip color zero */
2487 }
2488 else {
2489 if (colors < 16)
2490 color = (int)(coloriter & andcolor);
2491 else
2492 color = (int)(((coloriter - 1) % andcolor) + 1);
2493 }
2494 }
2495 if(debugflag != 470)
2496 if(color <= 0 && stdcalcmode == 'b' ) /* fix BTM bug */
2497 color = 1;
2498 (*plot) (col, row, color);
2499
2500 maxit = savemaxit;
2501 if ((kbdcount -= abs((int)realcoloriter)) <= 0)
2502 {
2503 if (check_key())
2504 return (-1);
2505 kbdcount = max_kbdcount;
2506 }
2507 return (color);
2508 }
2509 #undef green
2510 #undef yellow
2511
2512 #define cos45 sin45
2513 #define lcos45 lsin45
2514
2515 /**************** standardfractal doodad subroutines *********************/
2516 static void decomposition(void)
2517 {
2518 /* static double cos45 = 0.70710678118654750; */ /* cos 45 degrees */
2519 static double sin45 = 0.70710678118654750; /* sin 45 degrees */
2520 static double cos22_5 = 0.92387953251128670; /* cos 22.5 degrees */
2521 static double sin22_5 = 0.38268343236508980; /* sin 22.5 degrees */
2522 static double cos11_25 = 0.98078528040323040; /* cos 11.25 degrees */
2523 static double sin11_25 = 0.19509032201612820; /* sin 11.25 degrees */
2524 static double cos5_625 = 0.99518472667219690; /* cos 5.625 degrees */
2525 static double sin5_625 = 0.09801714032956060; /* sin 5.625 degrees */
2526 static double tan22_5 = 0.41421356237309500; /* tan 22.5 degrees */
2527 static double tan11_25 = 0.19891236737965800; /* tan 11.25 degrees */
2528 static double tan5_625 = 0.09849140335716425; /* tan 5.625 degrees */
2529 static double tan2_8125 = 0.04912684976946725; /* tan 2.8125 degrees */
2530 static double tan1_4063 = 0.02454862210892544; /* tan 1.4063 degrees */
2531 /* static long lcos45 ;*/ /* cos 45 degrees */
2532 static long lsin45 ; /* sin 45 degrees */
2533 static long lcos22_5 ; /* cos 22.5 degrees */
2534 static long lsin22_5 ; /* sin 22.5 degrees */
2535 static long lcos11_25 ; /* cos 11.25 degrees */
2536 static long lsin11_25 ; /* sin 11.25 degrees */
2537 static long lcos5_625 ; /* cos 5.625 degrees */
2538 static long lsin5_625 ; /* sin 5.625 degrees */
2539 static long ltan22_5 ; /* tan 22.5 degrees */
2540 static long ltan11_25 ; /* tan 11.25 degrees */
2541 static long ltan5_625 ; /* tan 5.625 degrees */
2542 static long ltan2_8125 ; /* tan 2.8125 degrees */
2543 static long ltan1_4063 ; /* tan 1.4063 degrees */
2544 static long reset_fudge = -1;
2545 int temp = 0;
2546 int save_temp = 0;
2547 int i;
2548 _LCMPLX lalt;
2549 _CMPLX alt;
2550 coloriter = 0;
2551 if (integerfractal) /* the only case */
2552 {
2553 if (reset_fudge != fudge)
2554 {
2555 reset_fudge = fudge;
2556 /* lcos45 = (long)(cos45 * fudge); */
2557 lsin45 = (long)(sin45 * fudge);
2558 lcos22_5 = (long)(cos22_5 * fudge);
2559 lsin22_5 = (long)(sin22_5 * fudge);
2560 lcos11_25 = (long)(cos11_25 * fudge);
2561 lsin11_25 = (long)(sin11_25 * fudge);
2562 lcos5_625 = (long)(cos5_625 * fudge);
2563 lsin5_625 = (long)(sin5_625 * fudge);
2564 ltan22_5 = (long)(tan22_5 * fudge);
2565 ltan11_25 = (long)(tan11_25 * fudge);
2566 ltan5_625 = (long)(tan5_625 * fudge);
2567 ltan2_8125 = (long)(tan2_8125 * fudge);
2568 ltan1_4063 = (long)(tan1_4063 * fudge);
2569 }
2570 if (lnew.y < 0)
2571 {
2572 temp = 2;
2573 lnew.y = -lnew.y;
2574 }
2575
2576 if (lnew.x < 0)
2577 {
2578 ++temp;
2579 lnew.x = -lnew.x;
2580 }
2581 if (decomp[0] == 2 && save_release >= 1827)
2582 {
2583 save_temp = temp;
2584 if(temp==2) save_temp = 3;
2585 if(temp==3) save_temp = 2;
2586 }
2587
2588 if (decomp[0] >= 8)
2589 {
2590 temp <<= 1;
2591 if (lnew.x < lnew.y)
2592 {
2593 ++temp;
2594 lalt.x = lnew.x; /* just */
2595 lnew.x = lnew.y; /* swap */
2596 lnew.y = lalt.x; /* them */
2597 }
2598
2599 if (decomp[0] >= 16)
2600 {
2601 temp <<= 1;
2602 if (multiply(lnew.x,ltan22_5,bitshift) < lnew.y)
2603 {
2604 ++temp;
2605 lalt = lnew;
2606 lnew.x = multiply(lalt.x,lcos45,bitshift) +
2607 multiply(lalt.y,lsin45,bitshift);
2608 lnew.y = multiply(lalt.x,lsin45,bitshift) -
2609 multiply(lalt.y,lcos45,bitshift);
2610 }
2611
2612 if (decomp[0] >= 32)
2613 {
2614 temp <<= 1;
2615 if (multiply(lnew.x,ltan11_25,bitshift) < lnew.y)
2616 {
2617 ++temp;
2618 lalt = lnew;
2619 lnew.x = multiply(lalt.x,lcos22_5,bitshift) +
2620 multiply(lalt.y,lsin22_5,bitshift);
2621 lnew.y = multiply(lalt.x,lsin22_5,bitshift) -
2622 multiply(lalt.y,lcos22_5,bitshift);
2623 }
2624
2625 if (decomp[0] >= 64)
2626 {
2627 temp <<= 1;
2628 if (multiply(lnew.x,ltan5_625,bitshift) < lnew.y)
2629 {
2630 ++temp;
2631 lalt = lnew;
2632 lnew.x = multiply(lalt.x,lcos11_25,bitshift) +
2633 multiply(lalt.y,lsin11_25,bitshift);
2634 lnew.y = multiply(lalt.x,lsin11_25,bitshift) -
2635 multiply(lalt.y,lcos11_25,bitshift);
2636 }
2637
2638 if (decomp[0] >= 128)
2639 {
2640 temp <<= 1;
2641 if (multiply(lnew.x,ltan2_8125,bitshift) < lnew.y)
2642 {
2643 ++temp;
2644 lalt = lnew;
2645 lnew.x = multiply(lalt.x,lcos5_625,bitshift) +
2646 multiply(lalt.y,lsin5_625,bitshift);
2647 lnew.y = multiply(lalt.x,lsin5_625,bitshift) -
2648 multiply(lalt.y,lcos5_625,bitshift);
2649 }
2650
2651 if (decomp[0] == 256)
2652 {
2653 temp <<= 1;
2654 if (multiply(lnew.x,ltan1_4063,bitshift) < lnew.y)
2655 if ((lnew.x*ltan1_4063 < lnew.y))
2656 ++temp;
2657 }
2658 }
2659 }
2660 }
2661 }
2662 }
2663 }
2664 else /* double case */
2665 {
2666 if (new.y < 0)
2667 {
2668 temp = 2;
2669 new.y = -new.y;
2670 }
2671 if (new.x < 0)
2672 {
2673 ++temp;
2674 new.x = -new.x;
2675 }
2676 if (decomp[0] == 2 && save_release >= 1827)
2677 {
2678 save_temp = temp;
2679 if(temp==2) save_temp = 3;
2680 if(temp==3) save_temp = 2;
2681 }
2682 if (decomp[0] >= 8)
2683 {
2684 temp <<= 1;
2685 if (new.x < new.y)
2686 {
2687 ++temp;
2688 alt.x = new.x; /* just */
2689 new.x = new.y; /* swap */
2690 new.y = alt.x; /* them */
2691 }
2692 if (decomp[0] >= 16)
2693 {
2694 temp <<= 1;
2695 if (new.x*tan22_5 < new.y)
2696 {
2697 ++temp;
2698 alt = new;
2699 new.x = alt.x*cos45 + alt.y*sin45;
2700 new.y = alt.x*sin45 - alt.y*cos45;
2701 }
2702
2703 if (decomp[0] >= 32)
2704 {
2705 temp <<= 1;
2706 if (new.x*tan11_25 < new.y)
2707 {
2708 ++temp;
2709 alt = new;
2710 new.x = alt.x*cos22_5 + alt.y*sin22_5;
2711 new.y = alt.x*sin22_5 - alt.y*cos22_5;
2712 }
2713
2714 if (decomp[0] >= 64)
2715 {
2716 temp <<= 1;
2717 if (new.x*tan5_625 < new.y)
2718 {
2719 ++temp;
2720 alt = new;
2721 new.x = alt.x*cos11_25 + alt.y*sin11_25;
2722 new.y = alt.x*sin11_25 - alt.y*cos11_25;
2723 }
2724
2725 if (decomp[0] >= 128)
2726 {
2727 temp <<= 1;
2728 if (new.x*tan2_8125 < new.y)
2729 {
2730 ++temp;
2731 alt = new;
2732 new.x = alt.x*cos5_625 + alt.y*sin5_625;
2733 new.y = alt.x*sin5_625 - alt.y*cos5_625;
2734 }
2735
2736 if (decomp[0] == 256)
2737 {
2738 temp <<= 1;
2739 if ((new.x*tan1_4063 < new.y))
2740 ++temp;
2741 }
2742 }
2743 }
2744 }
2745 }
2746 }
2747 }
2748 for (i = 1; temp > 0; ++i)
2749 {
2750 if (temp & 1)
2751 coloriter = (1 << i) - 1 - coloriter;
2752 temp >>= 1;
2753 }
2754 if (decomp[0] == 2 && save_release >= 1827) {
2755 if(save_temp & 2) coloriter = 1;
2756 else coloriter = 0;
2757 if (colors == 2)
2758 coloriter++;
2759 }
2760 else if (decomp[0] == 2 && save_release < 1827)
2761 coloriter &= 1;
2762 if (colors > decomp[0])
2763 coloriter++;
2764 }
2765
2766 /******************************************************************/
2767 /* Continuous potential calculation for Mandelbrot and Julia */
2768 /* Reference: Science of Fractal Images p. 190. */
2769 /* Special thanks to Mark Peterson for his "MtMand" program that */
2770 /* beautifully approximates plate 25 (same reference) and spurred */
2771 /* on the inclusion of similar capabilities in FRACTINT. */
2772 /* */
2773 /* The purpose of this function is to calculate a color value */
2774 /* for a fractal that varies continuously with the screen pixels */
2775 /* locations for better rendering in 3D. */
2776 /* */
2777 /* Here "magnitude" is the modulus of the orbit value at */
2778 /* "iterations". The potparms[] are user-entered paramters */
2779 /* controlling the level and slope of the continuous potential */
2780 /* surface. Returns color. - Tim Wegner 6/25/89 */
2781 /* */
2782 /* -- Change history -- */
2783 /* */
2784 /* 09/12/89 - added floatflag support and fixed float underflow */
2785 /* */
2786 /******************************************************************/
2787
2788 static int _fastcall potential(double mag, long iterations)
2789 {
2790 float f_mag,f_tmp,pot;
2791 double d_tmp;
2792 int i_pot;
2793 long l_pot;
2794
2795 if(iterations < maxit)
2796 {
2797 pot = (float)(l_pot = iterations+2);
2798 if(l_pot <= 0 || mag <= 1.0)
2799 pot = (float)0.0;
2800 else /* pot = log(mag) / pow(2.0, (double)pot); */
2801 {
2802 if(l_pot < 120 && !floatflag) /* empirically determined limit of fShift */
2803 {
2804 f_mag = (float)mag;
2805 fLog14(f_mag,f_tmp); /* this SHOULD be non-negative */
2806 fShift(f_tmp,(char)-l_pot,pot);
2807 }
2808 else
2809 {
2810 d_tmp = log(mag)/(double)pow(2.0,(double)pot);
2811 if(d_tmp > FLT_MIN) /* prevent float type underflow */
2812 pot = (float)d_tmp;
2813 else
2814 pot = (float)0.0;
2815 }
2816 }
2817 /* following transformation strictly for aesthetic reasons */
2818 /* meaning of parameters:
2819 potparam[0] -- zero potential level - highest color -
2820 potparam[1] -- slope multiplier -- higher is steeper
2821 potparam[2] -- rqlim value if changeable (bailout for modulus) */
2822
2823 if(pot > 0.0)
2824 {
2825 if(floatflag)
2826 pot = (float)sqrt((double)pot);
2827 else
2828 {
2829 fSqrt14(pot,f_tmp);
2830 pot = f_tmp;
2831 }
2832 pot = (float)(potparam[0] - pot*potparam[1] - 1.0);
2833 }
2834 else
2835 pot = (float)(potparam[0] - 1.0);
2836 if(pot < 1.0)
2837 pot = (float)1.0; /* avoid color 0 */
2838 }
2839 else if(inside >= 0)
2840 pot = inside;
2841 else /* inside < 0 implies inside=maxit, so use 1st pot param instead */
2842 pot = (float)potparam[0];
2843
2844 i_pot = (int)((l_pot = (long)(pot * 256)) >> 8);
2845 if(i_pot >= colors)
2846 {
2847 i_pot = colors - 1;
2848 l_pot = 255;
2849 }
2850
2851 if(pot16bit)
2852 {
2853 if (dotmode != 11) /* if putcolor won't be doing it for us */
2854 writedisk(col+sxoffs,row+syoffs,i_pot);
2855 writedisk(col+sxoffs,row+sydots+syoffs,(int)l_pot);
2856 }
2857
2858 return(i_pot);
2859 }
2860
2861
2862 /******************* boundary trace method ***************************
2863 Fractint's original btm was written by David Guenther. There were a few
2864 rare circumstances in which the original btm would not trace or fill
2865 correctly, even on Mandelbrot Sets. The code below was adapted from
2866 "Mandelbrot Sets by Wesley Loewer" (see calmanfp.asm) which was written
2867 before I was introduced to Fractint. It should be noted that without
2868 David Guenther's implimentation of a btm, I doubt that I would have been
2869 able to impliment my own code into Fractint. There are several things in
2870 the following code that are not original with me but came from David
2871 Guenther's code. I've noted these places with the initials DG.
2872
2873 Wesley Loewer 3/8/92
2874 *********************************************************************/
2875 #define bkcolor 0 /* I have some ideas for the future with this. -Wes */
2876 #define advance_match() coming_from = ((going_to = (going_to - 1) & 0x03) - 1) & 0x03
2877 #define advance_no_match() going_to = (going_to + 1) & 0x03
2878
2879 static
2880 int bound_trace_main(void)
2881 {
2882 enum direction coming_from;
2883 unsigned int match_found, continue_loop;
2884 int trail_color, fillcolor_used, last_fillcolor_used = -1;
2885 int max_putline_length;
2886 int right, left, length;
2887 static FCODE btm_cantbeused[]={"Boundary tracing cannot be used with "};
2888 if (inside == 0 || outside == 0)
2889 {
2890 static FCODE inside_outside[] = {"inside=0 or outside=0"};
2891 char msg[MSGLEN];
2892 far_strcpy(msg,btm_cantbeused);
2893 far_strcat(msg,inside_outside);
2894 stopmsg(0,msg);
2895 return(-1);
2896 }
2897 if (colors < 16)
2898 {
2899 char msg[MSGLEN];
2900 static FCODE lessthansixteen[] = {"< 16 colors"};
2901 far_strcpy(msg,btm_cantbeused);
2902 far_strcat(msg,lessthansixteen);
2903 stopmsg(0,msg);
2904 return(-1);
2905 }
2906
2907 got_status = 2;
2908 max_putline_length = 0; /* reset max_putline_length */
2909 for (currow = iystart; currow <= iystop; currow++)
2910 {
2911 reset_periodicity = 1; /* reset for a new row */
2912 color = bkcolor;
2913 for (curcol = ixstart; curcol <= ixstop; curcol++)
2914 {
2915 if (getcolor(curcol, currow) != bkcolor)
2916 continue;
2917
2918 trail_color = color;
2919 row = currow;
2920 col = curcol;
2921 if ((*calctype)()== -1) /* color, row, col are global */
2922 {
2923 if (showdot != bkcolor) /* remove showdot pixel */
2924 (*plot)(col,row,bkcolor);
2925 if (iystop != yystop) /* DG */
2926 iystop = yystop - (currow - yystart); /* allow for sym */
2927 add_worklist(xxstart,xxstop,curcol,currow,iystop,currow,0,worksym);
2928 return -1;
2929 }
2930 reset_periodicity = 0; /* normal periodicity checking */
2931
2932 /*
2933 This next line may cause a few more pixels to be calculated,
2934 but at the savings of quite a bit of overhead
2935 */
2936 if (color != trail_color) /* DG */
2937 continue;
2938
2939 /* sweep clockwise to trace outline */
2940 trail_row = currow;
2941 trail_col = curcol;
2942 trail_color = color;
2943 fillcolor_used = fillcolor > 0 ? fillcolor : trail_color;
2944 coming_from = West;
2945 going_to = East;
2946 match_found = 0;
2947 continue_loop = TRUE;
2948 do
2949 {
2950 step_col_row();
2951 if (row >= currow
2952 && col >= ixstart
2953 && col <= ixstop
2954 && row <= iystop)
2955 {
2956 /* the order of operations in this next line is critical */
2957 if ((color = getcolor(col, row)) == bkcolor && (*calctype)()== -1)
2958 /* color, row, col are global for (*calctype)() */
2959 {
2960 if (showdot != bkcolor) /* remove showdot pixel */
2961 (*plot)(col,row,bkcolor);
2962 if (iystop != yystop) /* DG */
2963 iystop = yystop - (currow - yystart); /* allow for sym */
2964 add_worklist(xxstart,xxstop,curcol,currow,iystop,currow,0,worksym);
2965 return -1;
2966 }
2967 else if (color == trail_color)
2968 {
2969 if (match_found < 4) /* to keep it from overflowing */
2970 match_found++;
2971 trail_row = row;
2972 trail_col = col;
2973 advance_match();
2974 }
2975 else
2976 {
2977 advance_no_match();
2978 continue_loop = going_to != coming_from || match_found;
2979 }
2980 }
2981 else
2982 {
2983 advance_no_match();
2984 continue_loop = going_to != coming_from || match_found;
2985 }
2986 } while (continue_loop && (col != curcol || row != currow));
2987
2988 if (match_found <= 3) /* DG */
2989 { /* no hole */
2990 color = bkcolor;
2991 reset_periodicity = 1;
2992 continue;
2993 }
2994
2995 /*
2996 Fill in region by looping around again, filling lines to the left
2997 whenever going_to is South or West
2998 */
2999 trail_row = currow;
3000 trail_col = curcol;
3001 coming_from = West;
3002 going_to = East;
3003 do
3004 {
3005 match_found = FALSE;
3006 do
3007 {
3008 step_col_row();
3009 if (row >= currow
3010 && col >= ixstart
3011 && col <= ixstop
3012 && row <= iystop
3013 && getcolor(col,row) == trail_color)
3014 /* getcolor() must be last */
3015 {
3016 if (going_to == South
3017 || (going_to == West && coming_from != East))
3018 { /* fill a row, but only once */
3019 right = col;
3020 while (--right >= ixstart && (color = getcolor(right,row)) == trail_color)
3021 ; /* do nothing */
3022 if (color == bkcolor) /* check last color */
3023 {
3024 left = right;
3025 while (getcolor(--left,row) == bkcolor)
3026 /* Should NOT be possible for left < ixstart */
3027 ; /* do nothing */
3028 left++; /* one pixel too far */
3029 if (right == left) /* only one hole */
3030 (*plot)(left,row,fillcolor_used);
3031 else
3032 { /* fill the line to the left */
3033 length=right-left+1;
3034 if (fillcolor_used != last_fillcolor_used || length > max_putline_length)
3035 { /* only reset dstack if necessary */
3036 memset(dstack,fillcolor_used,length);
3037 last_fillcolor_used = fillcolor_used;
3038 max_putline_length = length;
3039 }
3040 sym_fill_line(row, left, right, dstack);
3041 }
3042 } /* end of fill line */
3043
3044 #if 0 /* don't interupt with a check_key() during fill */
3056 #endif
3057 }
3058 trail_row = row;
3059 trail_col = col;
3060 advance_match();
3061 match_found = TRUE;
3062 }
3063 else
3064 advance_no_match();
3065 } while (!match_found && going_to != coming_from);
3066
3067 if (!match_found)
3068 { /* next one has to be a match */
3069 step_col_row();
3070 trail_row = row;
3071 trail_col = col;
3072 advance_match();
3073 }
3074 } while (trail_col != curcol || trail_row != currow);
3075 reset_periodicity = 1; /* reset after a trace/fill */
3076 color = bkcolor;
3077 }
3078 }
3079 return 0;
3080 }
3081
3082 /*******************************************************************/
3083 /* take one step in the direction of going_to */
3084 static void step_col_row()
3085 {
3086 switch (going_to)
3087 {
3088 case North:
3089 col = trail_col;
3090 row = trail_row - 1;
3091 break;
3092 case East:
3093 col = trail_col + 1;
3094 row = trail_row;
3095 break;
3096 case South:
3097 col = trail_col;
3098 row = trail_row + 1;
3099 break;
3100 case West:
3101 col = trail_col - 1;
3102 row = trail_row;
3103 break;
3104 }
3105 }
3106
3107 /******************* end of boundary trace method *******************/
3108
3109
3110 /************************ super solid guessing *****************************/
3111
3112 /*
3113 I, Timothy Wegner, invented this solidguessing idea and implemented it in
3114 more or less the overall framework you see here. I am adding this note
3115 now in a possibly vain attempt to secure my place in history, because
3116 Pieter Branderhorst has totally rewritten this routine, incorporating
3117 a *MUCH* more sophisticated algorithm. His revised code is not only
3118 faster, but is also more accurate. Harrumph!
3119 */
3120
3121 static int solidguess(void)
3122 {
3123 int i,x,y,xlim,ylim,blocksize;
3124 unsigned int *pfxp0,*pfxp1;
3125 unsigned int u;
3126
3127 guessplot=(plot!=putcolor && plot!=symplot2 && plot!=symplot2J);
3128 /* check if guessing at bottom & right edges is ok */
3129 bottom_guess = (plot == symplot2 || (plot == putcolor && iystop+1 == ydots));
3130 right_guess = (plot == symplot2J
3131 || ((plot == putcolor || plot == symplot2) && ixstop+1 == xdots));
3132
3133 /* there seems to be a bug in solid guessing at bottom and side */
3134 if(debugflag != 472)
3135 bottom_guess = right_guess = 0; /* TIW march 1995 */
3136
3137 i = maxblock = blocksize = ssg_blocksize();
3138 totpasses = 1;
3139 while ((i >>= 1) > 1) ++totpasses;
3140
3141 /* ensure window top and left are on required boundary, treat window
3142 as larger than it really is if necessary (this is the reason symplot
3143 routines must check for > xdots/ydots before plotting sym points) */
3144 ixstart &= -1 - (maxblock-1);
3145 iystart = yybegin;
3146 iystart &= -1 - (maxblock-1);
3147
3148 got_status = 1;
3149
3150 if (workpass == 0) /* otherwise first pass already done */
3151 {
3152 /* first pass, calc every blocksize**2 pixel, quarter result & paint it */
3153 curpass = 1;
3154 if (iystart <= yystart) /* first time for this window, init it */
3155 {
3156 currow = 0;
3157 memset(&tprefix[1][0][0],0,maxxblk*maxyblk*2); /* noskip flags off */
3158 reset_periodicity = 1;
3159 row=iystart;
3160 for(col=ixstart; col<=ixstop; col+=maxblock)
3161 { /* calc top row */
3162 if((*calctype)()== -1)
3163 {
3164 add_worklist(xxstart,xxstop,xxbegin,yystart,yystop,yybegin,0,worksym);
3165 goto exit_solidguess;
3166 }
3167 reset_periodicity = 0;
3168 }
3169 }
3170 else
3171 memset(&tprefix[1][0][0],-1,maxxblk*maxyblk*2); /* noskip flags on */
3172 for(y=iystart; y<=iystop; y+=blocksize)
3173 {
3174 currow = y;
3175 i = 0;
3176 if(y+blocksize<=iystop)
3177 { /* calc the row below */
3178 row=y+blocksize;
3179 reset_periodicity = 1;
3180 for(col=ixstart; col<=ixstop; col+=maxblock)
3181 {
3182 if((i=(*calctype)()) == -1)
3183 break;
3184 reset_periodicity = 0;
3185 }
3186 }
3187 reset_periodicity = 0;
3188 if (i == -1 || guessrow(1,y,blocksize) != 0) /* interrupted? */
3189 {
3190 if (y < yystart)
3191 y = yystart;
3192 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,y,0,worksym);
3193 goto exit_solidguess;
3194 }
3195 }
3196
3197 if (num_worklist) /* work list not empty, just do 1st pass */
3198 {
3199 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,yystart,1,worksym);
3200 goto exit_solidguess;
3201 }
3202 ++workpass;
3203 iystart = yystart & (-1 - (maxblock-1));
3204
3205 /* calculate skip flags for skippable blocks */
3206 xlim=(ixstop+maxblock)/maxblock+1;
3207 ylim=((iystop+maxblock)/maxblock+15)/16+1;
3208 if(right_guess==0) /* no right edge guessing, zap border */
3209 for(y=0;y<=ylim;++y)
3210 tprefix[1][y][xlim]= 0xffff;
3211 if(bottom_guess==0) /* no bottom edge guessing, zap border */
3212 {
3213 i=(iystop+maxblock)/maxblock+1;
3214 y=i/16+1;
3215 i=1<<(i&15);
3216 for(x=0;x<=xlim;++x)
3217 tprefix[1][y][x]|=i;
3218 }
3219 /* set each bit in tprefix[0] to OR of it & surrounding 8 in tprefix[1] */
3220 for(y=0;++y<ylim;)
3221 {
3222 pfxp0= (unsigned int *)&tprefix[0][y][0];
3223 pfxp1= (unsigned int *)&tprefix[1][y][0];
3224 for(x=0;++x<xlim;)
3225 {
3226 ++pfxp1;
3227 u= *(pfxp1-1)|*pfxp1|*(pfxp1+1);
3228 *(++pfxp0)=u|(u>>1)|(u<<1)
3229 |((*(pfxp1-(maxxblk+1))|*(pfxp1-maxxblk)|*(pfxp1-(maxxblk-1)))>>15)
3230 |((*(pfxp1+(maxxblk-1))|*(pfxp1+maxxblk)|*(pfxp1+(maxxblk+1)))<<15);
3231 }
3232 }
3233 }
3234 else /* first pass already done */
3235 memset(&tprefix[0][0][0],-1,maxxblk*maxyblk*2); /* noskip flags on */
3236 if(three_pass)
3237 goto exit_solidguess;
3238
3239 /* remaining pass(es), halve blocksize & quarter each blocksize**2 */
3240 i = workpass;
3241 while (--i > 0) /* allow for already done passes */
3242 blocksize = blocksize>>1;
3243 reset_periodicity = 0;
3244 while((blocksize=blocksize>>1)>=2)
3245 {
3246 if(stoppass > 0)
3247 if(workpass >= stoppass)
3248 goto exit_solidguess;
3249 curpass = workpass + 1;
3250 for(y=iystart; y<=iystop; y+=blocksize)
3251 {
3252 currow = y;
3253 if(guessrow(0,y,blocksize) != 0)
3254 {
3255 if (y < yystart)
3256 y = yystart;
3257 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,y,workpass,worksym);
3258 goto exit_solidguess;
3259 }
3260 }
3261 ++workpass;
3262 if (num_worklist /* work list not empty, do one pass at a time */
3263 && blocksize>2) /* if 2, we just did last pass */
3264 {
3265 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,yystart,workpass,worksym);
3266 goto exit_solidguess;
3267 }
3268 iystart = yystart & (-1 - (maxblock-1));
3269 }
3270
3271 exit_solidguess:
3272 return(0);
3273 }
3274
3275 #define calcadot(c,x,y) { col=x; row=y; if((c=(*calctype)())== -1) return -1; }
3276
3277 static int _fastcall guessrow(int firstpass,int y,int blocksize)
3278 {
3279 int x,i,j,color;
3280 int xplushalf,xplusblock;
3281 int ylessblock,ylesshalf,yplushalf,yplusblock;
3282 int c21,c31,c41; /* cxy is the color of pixel at (x,y) */
3283 int c12,c22,c32,c42; /* where c22 is the topleft corner of */
3284 int c13,c23,c33; /* the block being handled in current */
3285 int c24, c44; /* iteration */
3286 int guessed23,guessed32,guessed33,guessed12,guessed13;
3287 int prev11,fix21,fix31;
3288 unsigned int *pfxptr,pfxmask;
3289
3290 c44 = c41 = c42 = 0; /* just for warning */
3291
3292 halfblock=blocksize>>1;
3293 i=y/maxblock;
3294 pfxptr= (unsigned int *)&tprefix[firstpass][(i>>4)+1][ixstart/maxblock];
3295 pfxmask=1<<(i&15);
3296 ylesshalf=y-halfblock;
3297 ylessblock=y-blocksize; /* constants, for speed */
3298 yplushalf=y+halfblock;
3299 yplusblock=y+blocksize;
3300 prev11= -1;
3301 c24=c12=c13=c22=getcolor(ixstart,y);
3302 c31=c21=getcolor(ixstart,(y>0)?ylesshalf:0);
3303 if(yplusblock<=iystop)
3304 c24=getcolor(ixstart,yplusblock);
3305 else if(bottom_guess==0)
3306 c24= -1;
3307 guessed12=guessed13=0;
3308
3309 for(x=ixstart; x<=ixstop;) /* increment at end, or when doing continue */
3310 {
3311 if((x&(maxblock-1))==0) /* time for skip flag stuff */
3312 {
3313 ++pfxptr;
3314 if(firstpass==0 && (*pfxptr&pfxmask)==0) /* check for fast skip */
3315 {
3316 /* next useful in testing to make skips visible */
3317 /*
3318 if(halfblock==1)
3319 {
3320 (*plot)(x+1,y,0); (*plot)(x,y+1,0); (*plot)(x+1,y+1,0);
3321 }
3322 */
3323 x+=maxblock;
3324 prev11=c31=c21=c24=c12=c13=c22;
3325 guessed12=guessed13=0;
3326 continue;
3327 }
3328 }
3329
3330 if(firstpass) /* 1st pass, paint topleft corner */
3331 plotblock(0,x,y,c22);
3332 /* setup variables */
3333 xplusblock=(xplushalf=x+halfblock)+halfblock;
3334 if(xplushalf>ixstop)
3335 {
3336 if(right_guess==0)
3337 c31= -1;
3338 }
3339 else if(y>0)
3340 c31=getcolor(xplushalf,ylesshalf);
3341 if(xplusblock<=ixstop)
3342 {
3343 if(yplusblock<=iystop)
3344 c44=getcolor(xplusblock,yplusblock);
3345 c41=getcolor(xplusblock,(y>0)?ylesshalf:0);
3346 c42=getcolor(xplusblock,y);
3347 }
3348 else if(right_guess==0)
3349 c41=c42=c44= -1;
3350 if(yplusblock>iystop)
3351 c44=(bottom_guess)?c42:-1;
3352
3353 /* guess or calc the remaining 3 quarters of current block */
3354 guessed23=guessed32=guessed33=1;
3355 c23=c32=c33=c22;
3356 if(yplushalf>iystop)
3357 {
3358 if(bottom_guess==0)
3359 c23=c33= -1;
3360 guessed23=guessed33= -1;
3361 guessed13=0; /* fix for ydots not divisible by four bug TW 2/16/97 */
3362 }
3363 if(xplushalf>ixstop)
3364 {
3365 if(right_guess==0)
3366 c32=c33= -1;
3367 guessed32=guessed33= -1;
3368 }
3369 for(;;) /* go around till none of 23,32,33 change anymore */
3370 {
3371 if(guessed33>0
3372 && (c33!=c44 || c33!=c42 || c33!=c24 || c33!=c32 || c33!=c23))
3373 {
3374 calcadot(c33,xplushalf,yplushalf);
3375 guessed33=0;
3376 }
3377 if(guessed32>0
3378 && (c32!=c33 || c32!=c42 || c32!=c31 || c32!=c21
3379 || c32!=c41 || c32!=c23))
3380 {
3381 calcadot(c32,xplushalf,y);
3382 guessed32=0;
3383 continue;
3384 }
3385 if(guessed23>0
3386 && (c23!=c33 || c23!=c24 || c23!=c13 || c23!=c12 || c23!=c32))
3387 {
3388 calcadot(c23,x,yplushalf);
3389 guessed23=0;
3390 continue;
3391 }
3392 break;
3393 }
3394
3395 if(firstpass) /* note whether any of block's contents were calculated */
3396 if(guessed23==0 || guessed32==0 || guessed33==0)
3397 *pfxptr|=pfxmask;
3398
3399 if(halfblock>1) { /* not last pass, check if something to display */
3400 if(firstpass) /* display guessed corners, fill in block */
3401 {
3402 if(guessplot)
3403 {
3404 if(guessed23>0)
3405 (*plot)(x,yplushalf,c23);
3406 if(guessed32>0)
3407 (*plot)(xplushalf,y,c32);
3408 if(guessed33>0)
3409 (*plot)(xplushalf,yplushalf,c33);
3410 }
3411 plotblock(1,x,yplushalf,c23);
3412 plotblock(0,xplushalf,y,c32);
3413 plotblock(1,xplushalf,yplushalf,c33);
3414 }
3415 else /* repaint changed blocks */
3416 {
3417 if(c23!=c22)
3418 plotblock(-1,x,yplushalf,c23);
3419 if(c32!=c22)
3420 plotblock(-1,xplushalf,y,c32);
3421 if(c33!=c22)
3422 plotblock(-1,xplushalf,yplushalf,c33);
3423 }
3424 }
3425
3426 /* check if some calcs in this block mean earlier guesses need fixing */
3427 fix21=((c22!=c12 || c22!=c32)
3428 && c21==c22 && c21==c31 && c21==prev11
3429 && y>0
3430 && (x==ixstart || c21==getcolor(x-halfblock,ylessblock))
3431 && (xplushalf>ixstop || c21==getcolor(xplushalf,ylessblock))
3432 && c21==getcolor(x,ylessblock));
3433 fix31=(c22!=c32
3434 && c31==c22 && c31==c42 && c31==c21 && c31==c41
3435 && y>0 && xplushalf<=ixstop
3436 && c31==getcolor(xplushalf,ylessblock)
3437 && (xplusblock>ixstop || c31==getcolor(xplusblock,ylessblock))
3438 && c31==getcolor(x,ylessblock));
3439 prev11=c31; /* for next time around */
3440 if(fix21)
3441 {
3442 calcadot(c21,x,ylesshalf);
3443 if(halfblock>1 && c21!=c22)
3444 plotblock(-1,x,ylesshalf,c21);
3445 }
3446 if(fix31)
3447 {
3448 calcadot(c31,xplushalf,ylesshalf);
3449 if(halfblock>1 && c31!=c22)
3450 plotblock(-1,xplushalf,ylesshalf,c31);
3451 }
3452 if(c23!=c22)
3453 {
3454 if(guessed12)
3455 {
3456 calcadot(c12,x-halfblock,y);
3457 if(halfblock>1 && c12!=c22)
3458 plotblock(-1,x-halfblock,y,c12);
3459 }
3460 if(guessed13)
3461 {
3462 calcadot(c13,x-halfblock,yplushalf);
3463 if(halfblock>1 && c13!=c22)
3464 plotblock(-1,x-halfblock,yplushalf,c13);
3465 }
3466 }
3467 c22=c42;
3468 c24=c44;
3469 c13=c33;
3470 c31=c21=c41;
3471 c12=c32;
3472 guessed12=guessed32;
3473 guessed13=guessed33;
3474 x+=blocksize;
3475 } /* end x loop */
3476
3477 if(firstpass==0 || guessplot) return 0;
3478
3479 /* paint rows the fast way */
3480 for(i=0;i<halfblock;++i)
3481 {
3482 if((j=y+i)<=iystop)
3483 put_line(j,xxstart,ixstop,&dstack[xxstart]);
3484 if((j=y+i+halfblock)<=iystop)
3485 put_line(j,xxstart,ixstop,&dstack[xxstart+OLDMAXPIXELS]);
3486 if(keypressed()) return -1;
3487 }
3488 if(plot!=putcolor) /* symmetry, just vertical & origin the fast way */
3489 {
3490 if(plot==symplot2J) /* origin sym, reverse lines */
3491 for(i=(ixstop+xxstart+1)/2;--i>=xxstart;)
3492 {
3493 color=dstack[i];
3494 dstack[i]=dstack[j=ixstop-(i-xxstart)];
3495 dstack[j]=(BYTE)color;
3496 j+=OLDMAXPIXELS;
3497 color=dstack[i+OLDMAXPIXELS];
3498 dstack[i+OLDMAXPIXELS]=dstack[j];
3499 dstack[j]=(BYTE)color;
3500 }
3501 for(i=0;i<halfblock;++i)
3502 {
3503 if((j=yystop-(y+i-yystart))>iystop && j<ydots)
3504 put_line(j,xxstart,ixstop,&dstack[xxstart]);
3505 if((j=yystop-(y+i+halfblock-yystart))>iystop && j<ydots)
3506 put_line(j,xxstart,ixstop,&dstack[xxstart+OLDMAXPIXELS]);
3507 if(keypressed()) return -1;
3508 }
3509 }
3510 return 0;
3511 }
3512
3513 static void _fastcall plotblock(int buildrow,int x,int y,int color)
3514 {
3515 int i,xlim,ylim;
3516 if((xlim=x+halfblock)>ixstop)
3517 xlim=ixstop+1;
3518 if(buildrow>=0 && guessplot==0) /* save it for later put_line */
3519 {
3520 if(buildrow==0)
3521 for(i=x;i<xlim;++i)
3522 dstack[i]=(BYTE)color;
3523 else
3524 for(i=x;i<xlim;++i)
3525 dstack[i+OLDMAXPIXELS]=(BYTE)color;
3526 if (x>=xxstart) /* when x reduced for alignment, paint those dots too */
3527 return; /* the usual case */
3528 }
3529 /* paint it */
3530 if((ylim=y+halfblock)>iystop)
3531 {
3532 if(y>iystop)
3533 return;
3534 ylim=iystop+1;
3535 }
3536 for(i=x;++i<xlim;)
3537 (*plot)(i,y,color); /* skip 1st dot on 1st row */
3538 while(++y<ylim)
3539 for(i=x;i<xlim;++i)
3540 (*plot)(i,y,color);
3541 }
3542
3543
3544 /************************* symmetry plot setup ************************/
3545
3546 static int _fastcall xsym_split(int xaxis_row,int xaxis_between)
3547 {
3548 int i;
3549 if ((worksym&0x11) == 0x10) /* already decided not sym */
3550 return(1);
3551 if ((worksym&1) != 0) /* already decided on sym */
3552 iystop = (yystart+yystop)/2;
3553 else /* new window, decide */
3554 {
3555 worksym |= 0x10;
3556 if (xaxis_row <= yystart || xaxis_row >= yystop)
3557 return(1); /* axis not in window */
3558 i = xaxis_row + (xaxis_row - yystart);
3559 if (xaxis_between)
3560 ++i;
3561 if (i > yystop) /* split into 2 pieces, bottom has the symmetry */
3562 {
3563 if (num_worklist >= MAXCALCWORK-1) /* no room to split */
3564 return(1);
3565 iystop = xaxis_row - (yystop - xaxis_row);
3566 if (!xaxis_between)
3567 --iystop;
3568 add_worklist(xxstart,xxstop,xxstart,iystop+1,yystop,iystop+1,workpass,0);
3569 yystop = iystop;
3570 return(1); /* tell set_symmetry no sym for current window */
3571 }
3572 if (i < yystop) /* split into 2 pieces, top has the symmetry */
3573 {
3574 if (num_worklist >= MAXCALCWORK-1) /* no room to split */
3575 return(1);
3576 add_worklist(xxstart,xxstop,xxstart,i+1,yystop,i+1,workpass,0);
3577 yystop = i;
3578 }
3579 iystop = xaxis_row;
3580 worksym |= 1;
3581 }
3582 symmetry = 0;
3583 return(0); /* tell set_symmetry its a go */
3584 }
3585
3586 static int _fastcall ysym_split(int yaxis_col,int yaxis_between)
3587 {
3588 int i;
3589 if ((worksym&0x22) == 0x20) /* already decided not sym */
3590 return(1);
3591 if ((worksym&2) != 0) /* already decided on sym */
3592 ixstop = (xxstart+xxstop)/2;
3593 else /* new window, decide */
3594 {
3595 worksym |= 0x20;
3596 if (yaxis_col <= xxstart || yaxis_col >= xxstop)
3597 return(1); /* axis not in window */
3598 i = yaxis_col + (yaxis_col - xxstart);
3599 if (yaxis_between)
3600 ++i;
3601 if (i > xxstop) /* split into 2 pieces, right has the symmetry */
3602 {
3603 if (num_worklist >= MAXCALCWORK-1) /* no room to split */
3604 return(1);
3605 ixstop = yaxis_col - (xxstop - yaxis_col);
3606 if (!yaxis_between)
3607 --ixstop;
3608 add_worklist(ixstop+1,xxstop,ixstop+1,yystart,yystop,yystart,workpass,0);
3609 xxstop = ixstop;
3610 return(1); /* tell set_symmetry no sym for current window */
3611 }
3612 if (i < xxstop) /* split into 2 pieces, left has the symmetry */
3613 {
3614 if (num_worklist >= MAXCALCWORK-1) /* no room to split */
3615 return(1);
3616 add_worklist(i+1,xxstop,i+1,yystart,yystop,yystart,workpass,0);
3617 xxstop = i;
3618 }
3619 ixstop = yaxis_col;
3620 worksym |= 2;
3621 }
3622 symmetry = 0;
3623 return(0); /* tell set_symmetry its a go */
3624 }
3625
3626 #ifdef _MSC_VER
3627 #pragma optimize ("ea", off)
3628 #endif
3629
3630 #if (_MSC_VER >= 700)
3631 #pragma code_seg ("calcfra1_text") /* place following in an overlay */
3632 #endif
3633
3634 static void _fastcall setsymmetry(int sym, int uselist) /* set up proper symmetrical plot functions */
3635 {
3636 int i;
3637 int parmszero, parmsnoreal, parmsnoimag;
3638 int xaxis_row, yaxis_col; /* pixel number for origin */
3639 int xaxis_between=0, yaxis_between=0; /* if axis between 2 pixels, not on one */
3640 int xaxis_on_screen=0, yaxis_on_screen=0;
3641 double ftemp;
3642 bf_t bft1;
3643 int saved=0;
3644 symmetry = 1;
3645 if(stdcalcmode == 's' || stdcalcmode == 'o')
3646 return;
3647 if(sym == NOPLOT && forcesymmetry == 999)
3648 {
3649 plot = noplot;
3650 return;
3651 }
3652 /* NOTE: 16-bit potential disables symmetry */
3653 /* also any decomp= option and any inversion not about the origin */
3654 /* also any rotation other than 180deg and any off-axis stretch */
3655 if(bf_math)
3656 if(cmp_bf(bfxmin,bfx3rd) || cmp_bf(bfymin,bfy3rd))
3657 return;
3658 if ((potflag && pot16bit) || (invert && inversion[2] != 0.0)
3659 || decomp[0] != 0
3660 || xxmin!=xx3rd || yymin!=yy3rd)
3661 return;
3662 if(sym != XAXIS && sym != XAXIS_NOPARM && inversion[1] != 0.0 && forcesymmetry == 999)
3663 return;
3664 if(forcesymmetry < 999)
3665 sym = forcesymmetry;
3666 else if(forcesymmetry == 1000)
3667 forcesymmetry = sym; /* for backwards compatibility */
3668 else if(outside==REAL || outside==IMAG || outside==MULT || outside==SUM
3669 || outside==ATAN || bailoutest==Manr || outside==FMOD)
3670 return;
3671 else if(inside==FMODI || outside==TDIS)
3672 return;
3673 parmszero = (parm.x == 0.0 && parm.y == 0.0 && useinitorbit != 1);
3674 parmsnoreal = (parm.x == 0.0 && useinitorbit != 1);
3675 parmsnoimag = (parm.y == 0.0 && useinitorbit != 1);
3676 switch (fractype)
3677 { case LMANLAMFNFN: /* These need only P1 checked. */
3678 case FPMANLAMFNFN: /* P2 is used for a switch value */
3679 case LMANFNFN: /* These have NOPARM set in fractalp.c, */
3680 case FPMANFNFN: /* but it only applies to P1. */
3681 case FPMANDELZPOWER: /* or P2 is an exponent */
3682 case LMANDELZPOWER:
3683 case FPMANZTOZPLUSZPWR:
3684 case MARKSMANDEL:
3685 case MARKSMANDELFP:
3686 case MARKSJULIA:
3687 case MARKSJULIAFP:
3688 break;
3689 case FORMULA: /* Check P2, P3, P4 and P5 */
3690 case FFORMULA:
3691 parmszero = (parmszero && param[2] == 0.0 && param[3] == 0.0
3692 && param[4] == 0.0 && param[5] == 0.0
3693 && param[6] == 0.0 && param[7] == 0.0
3694 && param[8] == 0.0 && param[9] == 0.0);
3695 parmsnoreal = (parmsnoreal && param[2] == 0.0 && param[4] == 0.0
3696 && param[6] == 0.0 && param[8] == 0.0);
3697 parmsnoimag = (parmsnoimag && param[3] == 0.0 && param[5] == 0.0
3698 && param[7] == 0.0 && param[9] == 0.0);
3699 break;
3700 default: /* Check P2 for the rest */
3701 parmszero = (parmszero && parm2.x == 0.0 && parm2.y == 0.0);
3702 }
3703 xaxis_row = yaxis_col = -1;
3704 if(bf_math)
3705 {
3706 saved = save_stack();
3707 bft1 = alloc_stack(rbflength+2);
3708 xaxis_on_screen = (sign_bf(bfymin) != sign_bf(bfymax));
3709 yaxis_on_screen = (sign_bf(bfxmin) != sign_bf(bfxmax));
3710 }
3711 else
3712 {
3713 xaxis_on_screen = (sign(yymin) != sign(yymax));
3714 yaxis_on_screen = (sign(xxmin) != sign(xxmax));
3715 }
3716 if (xaxis_on_screen) /* axis is on screen */
3717 {
3718 if(bf_math)
3719 {
3720 /* ftemp = (0.0-yymax) / (yymin-yymax); */
3721 sub_bf(bft1,bfymin,bfymax);
3722 div_bf(bft1,bfymax,bft1);
3723 neg_a_bf(bft1);
3724 ftemp = (double)bftofloat(bft1);
3725 }
3726 else
3727 ftemp = (0.0-yymax) / (yymin-yymax);
3728 ftemp *= (ydots-1);
3729 ftemp += 0.25;
3730 xaxis_row = (int)ftemp;
3731 xaxis_between = (ftemp - xaxis_row >= 0.5);
3732 if (uselist == 0 && (!xaxis_between || (xaxis_row+1)*2 != ydots))
3733 xaxis_row = -1; /* can't split screen, so dead center or not at all */
3734 }
3735 if (yaxis_on_screen) /* axis is on screen */
3736 {
3737 if(bf_math)
3738 {
3739 /* ftemp = (0.0-xxmin) / (xxmax-xxmin); */
3740 sub_bf(bft1,bfxmax,bfxmin);
3741 div_bf(bft1,bfxmin,bft1);
3742 neg_a_bf(bft1);
3743 ftemp = (double)bftofloat(bft1);
3744 }
3745 else
3746 ftemp = (0.0-xxmin) / (xxmax-xxmin);
3747 ftemp *= (xdots-1);
3748 ftemp += 0.25;
3749 yaxis_col = (int)ftemp;
3750 yaxis_between = (ftemp - yaxis_col >= 0.5);
3751 if (uselist == 0 && (!yaxis_between || (yaxis_col+1)*2 != xdots))
3752 yaxis_col = -1; /* can't split screen, so dead center or not at all */
3753 }
3754 switch(sym) /* symmetry switch */
3755 {
3756 case XAXIS_NOREAL: /* X-axis Symmetry (no real param) */
3757 if (!parmsnoreal) break;
3758 goto xsym;
3759 case XAXIS_NOIMAG: /* X-axis Symmetry (no imag param) */
3760 if (!parmsnoimag) break;
3761 goto xsym;
3762 case XAXIS_NOPARM: /* X-axis Symmetry (no params)*/
3763 if (!parmszero)
3764 break;
3765 xsym:
3766 case XAXIS: /* X-axis Symmetry */
3767 if (xsym_split(xaxis_row,xaxis_between) == 0) {
3768 if(basin)
3769 plot = symplot2basin;
3770 else
3771 plot = symplot2;
3772 }
3773 break;
3774 case YAXIS_NOPARM: /* Y-axis Symmetry (No Parms)*/
3775 if (!parmszero)
3776 break;
3777 case YAXIS: /* Y-axis Symmetry */
3778 if (ysym_split(yaxis_col,yaxis_between) == 0)
3779 plot = symplot2Y;
3780 break;
3781 case XYAXIS_NOPARM: /* X-axis AND Y-axis Symmetry (no parms)*/
3782 if(!parmszero)
3783 break;
3784 case XYAXIS: /* X-axis AND Y-axis Symmetry */
3785 xsym_split(xaxis_row,xaxis_between);
3786 ysym_split(yaxis_col,yaxis_between);
3787 switch (worksym & 3)
3788 {
3789 case 1: /* just xaxis symmetry */
3790 if(basin)
3791 plot = symplot2basin;
3792 else
3793 plot = symplot2 ;
3794 break;
3795 case 2: /* just yaxis symmetry */
3796 if (basin) /* got no routine for this case */
3797 {
3798 ixstop = xxstop; /* fix what split should not have done */
3799 symmetry = 1;
3800 }
3801 else
3802 plot = symplot2Y;
3803 break;
3804 case 3: /* both axes */
3805 if(basin)
3806 plot = symplot4basin;
3807 else
3808 plot = symplot4 ;
3809 }
3810 break;
3811 case ORIGIN_NOPARM: /* Origin Symmetry (no parms)*/
3812 if (!parmszero)
3813 break;
3814 case ORIGIN: /* Origin Symmetry */
3815 originsym:
3816 if ( xsym_split(xaxis_row,xaxis_between) == 0
3817 && ysym_split(yaxis_col,yaxis_between) == 0)
3818 {
3819 plot = symplot2J;
3820 ixstop = xxstop; /* didn't want this changed */
3821 }
3822 else
3823 {
3824 iystop = yystop; /* in case first split worked */
3825 symmetry = 1;
3826 worksym = 0x30; /* let it recombine with others like it */
3827 }
3828 break;
3829 case PI_SYM_NOPARM:
3830 if (!parmszero)
3831 break;
3832 case PI_SYM: /* PI symmetry */
3833 if(bf_math)
3834 {
3835 if((double)bftofloat(abs_a_bf(sub_bf(bft1,bfxmax,bfxmin))) < PI/4)
3836 break; /* no point in pi symmetry if values too close */
3837 }
3838 else
3839 {
3840 if(fabs(xxmax - xxmin) < PI/4)
3841 break; /* no point in pi symmetry if values too close */
3842 }
3843 if(invert && forcesymmetry == 999)
3844 goto originsym;
3845 plot = symPIplot ;
3846 symmetry = 0;
3847 if ( xsym_split(xaxis_row,xaxis_between) == 0
3848 && ysym_split(yaxis_col,yaxis_between) == 0)
3849 if(parm.y == 0.0)
3850 plot = symPIplot4J; /* both axes */
3851 else
3852 plot = symPIplot2J; /* origin */
3853 else
3854 {
3855 iystop = yystop; /* in case first split worked */
3856 worksym = 0x30; /* don't mark pisym as ysym, just do it unmarked */
3857 }
3858 if(bf_math)
3859 {
3860 sub_bf(bft1,bfxmax,bfxmin);
3861 abs_a_bf(bft1);
3862 pixelpi = (int)((PI/(double)bftofloat(bft1)*xdots)); /* PI in pixels */
3863 }
3864 else
3865 pixelpi = (int)((PI/fabs(xxmax-xxmin))*xdots); /* PI in pixels */
3866
3867 if ( (ixstop = xxstart+pixelpi-1 ) > xxstop)
3868 ixstop = xxstop;
3869 if (plot == symPIplot4J && ixstop > (i = (xxstart+xxstop)/2))
3870 ixstop = i;
3871 break;
3872 default: /* no symmetry */
3873 break;
3874 }
3875 if(bf_math)
3876 restore_stack(saved);
3877 }
3878
3879 #ifdef _MSC_VER
3880 #pragma optimize ("ea", on)
3881 #endif
3882
3883 #if (_MSC_VER >= 700)
3884 #pragma code_seg () /* back to normal segment */
3885 #endif
3886
3887 /**************** tesseral method by CJLT begins here*********************/
3888 /* reworked by PB for speed and resumeability */
3889
3890 struct tess { /* one of these per box to be done gets stacked */
3891 int x1,x2,y1,y2; /* left/right top/bottom x/y coords */
3892 int top,bot,lft,rgt; /* edge colors, -1 mixed, -2 unknown */
3893 };
3894
3895 static int tesseral(void)
3896 {
3897 struct tess *tp;
3898
3899 guessplot = (plot != putcolor && plot != symplot2);
3900 tp = (struct tess *)&dstack[0];
3901 tp->x1 = ixstart; /* set up initial box */
3902 tp->x2 = ixstop;
3903 tp->y1 = iystart;
3904 tp->y2 = iystop;
3905
3906 if (workpass == 0) { /* not resuming */
3907 tp->top = tessrow(ixstart,ixstop,iystart); /* Do top row */
3908 tp->bot = tessrow(ixstart,ixstop,iystop); /* Do bottom row */
3909 tp->lft = tesscol(ixstart,iystart+1,iystop-1); /* Do left column */
3910 tp->rgt = tesscol(ixstop,iystart+1,iystop-1); /* Do right column */
3911 if (check_key()) { /* interrupt before we got properly rolling */
3912 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,yystart,0,worksym);
3913 return(-1);
3914 }
3915 }
3916
3917 else { /* resuming, rebuild work stack */
3918 int i,mid,curx,cury,xsize,ysize;
3919 struct tess *tp2;
3920 tp->top = tp->bot = tp->lft = tp->rgt = -2;
3921 cury = yybegin & 0xfff;
3922 ysize = 1;
3923 i = (unsigned)yybegin >> 12;
3924 while (--i >= 0) ysize <<= 1;
3925 curx = workpass & 0xfff;
3926 xsize = 1;
3927 i = (unsigned)workpass >> 12;
3928 while (--i >= 0) xsize <<= 1;
3929 for(;;) {
3930 tp2 = tp;
3931 if (tp->x2 - tp->x1 > tp->y2 - tp->y1) { /* next divide down middle */
3932 if (tp->x1 == curx && (tp->x2 - tp->x1 - 2) < xsize)
3933 break;
3934 mid = (tp->x1 + tp->x2) >> 1; /* Find mid point */
3935 if (mid > curx) { /* stack right part */
3936 memcpy(++tp,tp2,sizeof(*tp));
3937 tp->x2 = mid;
3938 }
3939 tp2->x1 = mid;
3940 }
3941 else { /* next divide across */
3942 if (tp->y1 == cury && (tp->y2 - tp->y1 - 2) < ysize) break;
3943 mid = (tp->y1 + tp->y2) >> 1; /* Find mid point */
3944 if (mid > cury) { /* stack bottom part */
3945 memcpy(++tp,tp2,sizeof(*tp));
3946 tp->y2 = mid;
3947 }
3948 tp2->y1 = mid;
3949 }
3950 }
3951 }
3952
3953 got_status = 4; /* for tab_display */
3954
3955 while (tp >= (struct tess *)&dstack[0]) { /* do next box */
3956 curcol = tp->x1; /* for tab_display */
3957 currow = tp->y1;
3958
3959 if (tp->top == -1 || tp->bot == -1 || tp->lft == -1 || tp->rgt == -1)
3960 goto tess_split;
3961 /* for any edge whose color is unknown, set it */
3962 if (tp->top == -2)
3963 tp->top = tesschkrow(tp->x1,tp->x2,tp->y1);
3964 if (tp->top == -1)
3965 goto tess_split;
3966 if (tp->bot == -2)
3967 tp->bot = tesschkrow(tp->x1,tp->x2,tp->y2);
3968 if (tp->bot != tp->top)
3969 goto tess_split;
3970 if (tp->lft == -2)
3971 tp->lft = tesschkcol(tp->x1,tp->y1,tp->y2);
3972 if (tp->lft != tp->top)
3973 goto tess_split;
3974 if (tp->rgt == -2)
3975 tp->rgt = tesschkcol(tp->x2,tp->y1,tp->y2);
3976 if (tp->rgt != tp->top)
3977 goto tess_split;
3978
3979 {
3980 int mid,midcolor;
3981 if (tp->x2 - tp->x1 > tp->y2 - tp->y1) { /* divide down the middle */
3982 mid = (tp->x1 + tp->x2) >> 1; /* Find mid point */
3983 midcolor = tesscol(mid, tp->y1+1, tp->y2-1); /* Do mid column */
3984 if (midcolor != tp->top) goto tess_split;
3985 }
3986 else { /* divide across the middle */
3987 mid = (tp->y1 + tp->y2) >> 1; /* Find mid point */
3988 midcolor = tessrow(tp->x1+1, tp->x2-1, mid); /* Do mid row */
3989 if (midcolor != tp->top) goto tess_split;
3990 }
3991 }
3992
3993 { /* all 4 edges are the same color, fill in */
3994 int i,j;
3995 i = 0;
3996 if(fillcolor != 0)
3997 {
3998 if(fillcolor > 0)
3999 tp->top = fillcolor %colors;
4000 if (guessplot || (j = tp->x2 - tp->x1 - 1) < 2) { /* paint dots */
4001 for (col = tp->x1 + 1; col < tp->x2; col++)
4002 for (row = tp->y1 + 1; row < tp->y2; row++) {
4003 (*plot)(col,row,tp->top);
4004 if (++i > 500) {
4005 if (check_key()) goto tess_end;
4006 i = 0;
4007 }
4008 }
4009 }
4010 else { /* use put_line for speed */
4011 memset(&dstack[OLDMAXPIXELS],tp->top,j);
4012 for (row = tp->y1 + 1; row < tp->y2; row++) {
4013 put_line(row,tp->x1+1,tp->x2-1,&dstack[OLDMAXPIXELS]);
4014 if (plot != putcolor) /* symmetry */
4015 if ((j = yystop-(row-yystart)) > iystop && j < ydots)
4016 put_line(j,tp->x1+1,tp->x2-1,&dstack[OLDMAXPIXELS]);
4017 if (++i > 25) {
4018 if (check_key()) goto tess_end;
4019 i = 0;
4020 }
4021 }
4022 }
4023 }
4024 --tp;
4025 }
4026 continue;
4027
4028 tess_split:
4029 { /* box not surrounded by same color, sub-divide */
4030 int mid,midcolor;
4031 struct tess *tp2;
4032 if (tp->x2 - tp->x1 > tp->y2 - tp->y1) { /* divide down the middle */
4033 mid = (tp->x1 + tp->x2) >> 1; /* Find mid point */
4034 midcolor = tesscol(mid, tp->y1+1, tp->y2-1); /* Do mid column */
4035 if (midcolor == -3) goto tess_end;
4036 if (tp->x2 - mid > 1) { /* right part >= 1 column */
4037 if (tp->top == -1) tp->top = -2;
4038 if (tp->bot == -1) tp->bot = -2;
4039 tp2 = tp;
4040 if (mid - tp->x1 > 1) { /* left part >= 1 col, stack right */
4041 memcpy(++tp,tp2,sizeof(*tp));
4042 tp->x2 = mid;
4043 tp->rgt = midcolor;
4044 }
4045 tp2->x1 = mid;
4046 tp2->lft = midcolor;
4047 }
4048 else
4049 --tp;
4050 }
4051 else { /* divide across the middle */
4052 mid = (tp->y1 + tp->y2) >> 1; /* Find mid point */
4053 midcolor = tessrow(tp->x1+1, tp->x2-1, mid); /* Do mid row */
4054 if (midcolor == -3) goto tess_end;
4055 if (tp->y2 - mid > 1) { /* bottom part >= 1 column */
4056 if (tp->lft == -1) tp->lft = -2;
4057 if (tp->rgt == -1) tp->rgt = -2;
4058 tp2 = tp;
4059 if (mid - tp->y1 > 1) { /* top also >= 1 col, stack bottom */
4060 memcpy(++tp,tp2,sizeof(*tp));
4061 tp->y2 = mid;
4062 tp->bot = midcolor;
4063 }
4064 tp2->y1 = mid;
4065 tp2->top = midcolor;
4066 }
4067 else
4068 --tp;
4069 }
4070 }
4071
4072 }
4073
4074 tess_end:
4075 if (tp >= (struct tess *)&dstack[0]) { /* didn't complete */
4076 int i,xsize,ysize;
4077 xsize = ysize = 1;
4078 i = 2;
4079 while (tp->x2 - tp->x1 - 2 >= i) {
4080 i <<= 1;
4081 ++xsize;
4082 }
4083 i = 2;
4084 while (tp->y2 - tp->y1 - 2 >= i) {
4085 i <<= 1;
4086 ++ysize;
4087 }
4088 add_worklist(xxstart,xxstop,xxstart,yystart,yystop,
4089 (ysize<<12)+tp->y1,(xsize<<12)+tp->x1,worksym);
4090 return(-1);
4091 }
4092 return(0);
4093
4094 } /* tesseral */
4095
4096 static int _fastcall tesschkcol(int x,int y1,int y2)
4097 {
4098 int i;
4099 i = getcolor(x,++y1);
4100 while (--y2 > y1)
4101 if (getcolor(x,y2) != i) return -1;
4102 return i;
4103 }
4104
4105 static int _fastcall tesschkrow(int x1,int x2,int y)
4106 {
4107 int i;
4108 i = getcolor(x1,y);
4109 while (x2 > x1) {
4110 if (getcolor(x2,y) != i) return -1;
4111 --x2;
4112 }
4113 return i;
4114 }
4115
4116 static int _fastcall tesscol(int x,int y1,int y2)
4117 {
4118 int colcolor,i;
4119 col = x;
4120 row = y1;
4121 reset_periodicity = 1;
4122 colcolor = (*calctype)();
4123 reset_periodicity = 0;
4124 while (++row <= y2) { /* generate the column */
4125 if ((i = (*calctype)()) < 0) return -3;
4126 if (i != colcolor) colcolor = -1;
4127 }
4128 return colcolor;
4129 }
4130
4131 static int _fastcall tessrow(int x1,int x2,int y)
4132 {
4133 int rowcolor,i;
4134 row = y;
4135 col = x1;
4136 reset_periodicity = 1;
4137 rowcolor = (*calctype)();
4138 reset_periodicity = 0;
4139 while (++col <= x2) { /* generate the row */
4140 if ((i = (*calctype)()) < 0) return -3;
4141 if (i != rowcolor) rowcolor = -1;
4142 }
4143 return rowcolor;
4144 }
4145
4146 /* added for testing autologmap() */ /* CAE 9211 fixed missing comment */
4147 /* insert at end of CALCFRAC.C */
4148
4149 static long autologmap(void) /*RB*/
4150 { /* calculate round screen edges to avoid wasted colours in logmap */
4151 long mincolour;
4152 int lag;
4153 int xstop = xdots - 1; /* don't use symetry */
4154 int ystop = ydots - 1; /* don't use symetry */
4155 long old_maxit;
4156 mincolour=LONG_MAX;
4157 row=0;
4158 reset_periodicity = 0;
4159 old_maxit = maxit;
4160 for (col=0;col<xstop;col++) /* top row */
4161 {
4162 color=(*calctype)();
4163 if (color == -1) goto ack; /* key pressed, bailout */
4164 if (realcoloriter < mincolour) {
4165 mincolour=realcoloriter ;
4166 maxit = max(2,mincolour); /*speedup for when edges overlap lakes */
4167 }
4168 if ( col >=32 ) (*plot)(col-32,row,0);
4169 } /* these lines tidy up for BTM etc */
4170 for (lag=32;lag>0;lag--) (*plot)(col-lag,row,0);
4171
4172 col=xstop;
4173 for (row=0;row<ystop;row++) /* right side */
4174 {
4175 color=(*calctype)();
4176 if (color == -1) goto ack; /* key pressed, bailout */
4177 if (realcoloriter < mincolour) {
4178 mincolour=realcoloriter ;
4179 maxit = max(2,mincolour); /*speedup for when edges overlap lakes */
4180 }
4181 if ( row >=32 ) (*plot)(col,row-32,0);
4182 }
4183 for (lag=32;lag>0;lag--) (*plot)(col,row-lag,0);
4184
4185 col=0;
4186 for (row=0;row<ystop;row++) /* left side */
4187 {
4188 color=(*calctype)();
4189 if (color == -1) goto ack; /* key pressed, bailout */
4190 if (realcoloriter < mincolour) {
4191 mincolour=realcoloriter ;
4192 maxit = max(2,mincolour); /*speedup for when edges overlap lakes */
4193 }
4194 if ( row >=32 ) (*plot)(col,row-32,0);
4195 }
4196 for (lag=32;lag>0;lag--) (*plot)(col,row-lag,0);
4197
4198 row=ystop ;
4199 for (col=0;col<xstop;col++) /* bottom row */
4200 {
4201 color=(*calctype)();
4202 if (color == -1) goto ack; /* key pressed, bailout */
4203 if (realcoloriter < mincolour) {
4204 mincolour=realcoloriter ;
4205 maxit = max(2,mincolour); /*speedup for when edges overlap lakes */
4206 }
4207 if ( col >=32 ) (*plot)(col-32,row,0);
4208 }
4209 for (lag=32;lag>0;lag--) (*plot)(col-lag,row,0);
4210
4211 ack: /* bailout here if key is pressed */
4212 if (mincolour==2) resuming=1; /* insure autologmap not called again */
4213 maxit = old_maxit;
4214
4215 return ( mincolour );
4216 }
4217
4218 /* Symmetry plot for period PI */
4219 void _fastcall symPIplot(int x, int y, int color)
4220 {
4221 while(x <= xxstop)
4222 {
4223 putcolor(x, y, color) ;
4224 x += pixelpi;
4225 }
4226 }
4227 /* Symmetry plot for period PI plus Origin Symmetry */
4228 void _fastcall symPIplot2J(int x, int y, int color)
4229 {
4230 int i,j;
4231 while(x <= xxstop)
4232 {
4233 putcolor(x, y, color) ;
4234 if ((i=yystop-(y-yystart)) > iystop && i < ydots
4235 && (j=xxstop-(x-xxstart)) < xdots)
4236 putcolor(j, i, color) ;
4237 x += pixelpi;
4238 }
4239 }
4240 /* Symmetry plot for period PI plus Both Axis Symmetry */
4241 void _fastcall symPIplot4J(int x, int y, int color)
4242 {
4243 int i,j;
4244 while(x <= (xxstart+xxstop)/2)
4245 {
4246 j = xxstop-(x-xxstart);
4247 putcolor( x , y , color) ;
4248 if (j < xdots)
4249 putcolor(j , y , color) ;
4250 if ((i=yystop-(y-yystart)) > iystop && i < ydots)
4251 {
4252 putcolor(x , i , color) ;
4253 if (j < xdots)
4254 putcolor(j , i , color) ;
4255 }
4256 x += pixelpi;
4257 }
4258 }
4259
4260 /* Symmetry plot for X Axis Symmetry */
4261 void _fastcall symplot2(int x, int y, int color)
4262 {
4263 int i;
4264 putcolor(x, y, color) ;
4265 if ((i=yystop-(y-yystart)) > iystop && i < ydots)
4266 putcolor(x, i, color) ;
4267 }
4268
4269 /* Symmetry plot for Y Axis Symmetry */
4270 void _fastcall symplot2Y(int x, int y, int color)
4271 {
4272 int i;
4273 putcolor(x, y, color) ;
4274 if ((i=xxstop-(x-xxstart)) < xdots)
4275 putcolor(i, y, color) ;
4276 }
4277
4278 /* Symmetry plot for Origin Symmetry */
4279 void _fastcall symplot2J(int x, int y, int color)
4280 {
4281 int i,j;
4282 putcolor(x, y, color) ;
4283 if ((i=yystop-(y-yystart)) > iystop && i < ydots
4284 && (j=xxstop-(x-xxstart)) < xdots)
4285 putcolor(j, i, color) ;
4286 }
4287
4288 /* Symmetry plot for Both Axis Symmetry */
4289 void _fastcall symplot4(int x, int y, int color)
4290 {
4291 int i,j;
4292 j = xxstop-(x-xxstart);
4293 putcolor( x , y, color) ;
4294 if (j < xdots)
4295 putcolor( j , y, color) ;
4296 if ((i=yystop-(y-yystart)) > iystop && i < ydots)
4297 {
4298 putcolor( x , i, color) ;
4299 if (j < xdots)
4300 putcolor(j , i, color) ;
4301 }
4302 }
4303
4304 /* Symmetry plot for X Axis Symmetry - Striped Newtbasin version */
4305 void _fastcall symplot2basin(int x, int y, int color)
4306 {
4307 int i,stripe;
4308 putcolor(x, y, color) ;
4309 if(basin==2 && color > 8)
4310 stripe=8;
4311 else
4312 stripe = 0;
4313 if ((i=yystop-(y-yystart)) > iystop && i < ydots)
4314 {
4315 color -= stripe; /* reconstruct unstriped color */
4316 color = (degree+1-color)%degree+1; /* symmetrical color */
4317 color += stripe; /* add stripe */
4318 putcolor(x, i,color) ;
4319 }
4320 }
4321
4322 /* Symmetry plot for Both Axis Symmetry - Newtbasin version */
4323 void _fastcall symplot4basin(int x, int y, int color)
4324 {
4325 int i,j,color1,stripe;
4326 if(color == 0) /* assumed to be "inside" color */
4327 {
4328 symplot4(x, y, color);
4329 return;
4330 }
4331 if(basin==2 && color > 8)
4332 stripe = 8;
4333 else
4334 stripe = 0;
4335 color -= stripe; /* reconstruct unstriped color */
4336 color1 = degree/2+degree+2 - color;
4337 if(color < degree/2+2)
4338 color1 = degree/2+2 - color;
4339 else
4340 color1 = degree/2+degree+2 - color;
4341 j = xxstop-(x-xxstart);
4342 putcolor( x, y, color+stripe) ;
4343 if (j < xdots)
4344 putcolor( j, y, color1+stripe) ;
4345 if ((i=yystop-(y-yystart)) > iystop && i < ydots)
4346 {
4347 putcolor( x, i, stripe + (degree+1 - color)%degree+1) ;
4348 if (j < xdots)
4349 putcolor(j, i, stripe + (degree+1 - color1)%degree+1) ;
4350 }
4351 }
4352
4353 static void _fastcall puttruecolor_disk(int x, int y, int color)
4354 {
4355 putcolor_a(x,y,color);
4356
4357 targa_color(x, y, color);
4358 /* writedisk(x*=3,y,(BYTE)((realcoloriter ) & 0xff)); */ /* blue */
4359 /* writedisk(++x, y,(BYTE)((realcoloriter >> 8 ) & 0xff)); */ /* green */
4360 /* writedisk(x+1, y,(BYTE)((realcoloriter >> 16) & 0xff)); */ /* red */
4361
4362 }
4363
4364 /* Do nothing plot!!! */
4365 #ifdef __CLINT__
4367 #endif
4368
4369 void _fastcall noplot(int x,int y,int color)
4370 {
4371 x = y = color = 0; /* just for warning */
4372 }
4373