fracsubr.c

frasr200的win 版本源码(18.21),使用make文件,使用的vc版本较低,在我的环境下编译有问题! 很不错的分形程序代码!

   ddelmin = fabs(delxx);
   if (fabs(delxx2) > ddelmin)
      ddelmin = fabs(delxx2);
   if (fabs(delyy) > fabs(delyy2)) {
      if (fabs(delyy) < ddelmin)
	 ddelmin = fabs(delyy);
      }
   else
      if (fabs(delyy2) < ddelmin)
	 ddelmin = fabs(delyy2);
   delmin = fudgetolong(ddelmin);

   /* calculate factors which plot real values to screen co-ords */
   /* calcfrac.c plot_orbit routines have comments about this	 */
   ftemp = (0.0-delyy2) * delxx2 * dxsize * dysize
	   - (xxmax-xx3rd) * (yy3rd-yymax);
   plotmx1 = delxx2 * dxsize * dysize / ftemp;
   plotmx2 = (yy3rd-yymax) * dxsize / ftemp;
   plotmy1 = (0.0-delyy2) * dxsize * dysize / ftemp;
   plotmy2 = (xxmax-xx3rd) * dysize / ftemp;

}

static long _fastcall fudgetolong(double d)
{
   if ((d *= fudge) > 0) d += 0.5;
   else 		 d -= 0.5;
   return (long)d;
}

static double _fastcall fudgetodouble(long l)
{
   char buf[30];
   double d;
   sprintf(buf,"%.9g",(double)l / fudge);
#ifndef XFRACT
   sscanf(buf,"%lg",&d);
#else
   sscanf(buf,"%lf",&d);
#endif
   return d;
}


void adjust_corner()
{
   /* make edges very near vert/horiz exact, to ditch rounding errs and */
   /* to avoid problems when delta per axis makes too large a ratio	*/
   double ftemp,ftemp2;
   if( (ftemp=fabs(xx3rd-xxmin)) < (ftemp2=fabs(xxmax-xx3rd)) ) {
      if (ftemp*10000 < ftemp2 && yy3rd != yymax)
	 xx3rd = xxmin;
      }
   else if (ftemp2*10000 < ftemp && yy3rd != yymin)
      xx3rd = xxmax;
   if( (ftemp=fabs(yy3rd-yymin)) < (ftemp2=fabs(yymax-yy3rd)) ) {
      if (ftemp*10000 < ftemp2 && xx3rd != xxmax)
	 yy3rd = yymin;
      }
   else if (ftemp2*10000 < ftemp && xx3rd != xxmin)
      yy3rd = yymax;
}

static void _fastcall adjust_to_limits(double expand)
{  double cornerx[4],cornery[4];
   double lowx,highx,lowy,highy,limit,ftemp;
   double centerx,centery,adjx,adjy;
   int i;
   limit = 32767.99;
   if (bitshift >= 24) limit = 31.99;
   if (bitshift >= 29) limit = 3.99;
   centerx = (xxmin+xxmax)/2;
   centery = (yymin+yymax)/2;
   if (xxmin == centerx) { /* ohoh, infinitely thin, fix it */
      smallest_add(&xxmax);
      xxmin -= xxmax-centerx;
      }
   if (yymin == centery) {
      smallest_add(&yymax);
      yymin -= yymax-centery;
      }
   if (xx3rd == centerx)
      smallest_add(&xx3rd);
   if (yy3rd == centery)
      smallest_add(&yy3rd);
   /* setup array for easier manipulation */
   cornerx[0] = xxmin; cornerx[1] = xxmax;
   cornerx[2] = xx3rd; cornerx[3] = xxmin+(xxmax-xx3rd);
   cornery[0] = yymax; cornery[1] = yymin;
   cornery[2] = yy3rd; cornery[3] = yymin+(yymax-yy3rd);
   /* if caller wants image size adjusted, do that first */
   if (expand != 1.0)
      for (i=0; i<4; ++i) {
	 cornerx[i] = centerx + (cornerx[i]-centerx)*expand;
	 cornery[i] = centery + (cornery[i]-centery)*expand;
      }
   /* get min/max x/y values */
   lowx = highx = cornerx[0];
   lowy = highy = cornery[0];
   for (i=1; i<4; ++i) {
      if (cornerx[i] < lowx)  lowx  = cornerx[i];
      if (cornerx[i] > highx) highx = cornerx[i];
      if (cornery[i] < lowy)  lowy  = cornery[i];
      if (cornery[i] > highy) highy = cornery[i];
      }
   /* if image is too large, downsize it maintaining center */
   ftemp = highx-lowx;
   if (highy-lowy > ftemp) ftemp = highy-lowy;
   if ((ftemp = limit*2/ftemp) < 1.0)
      for (i=0; i<4; ++i) {
	 cornerx[i] = centerx + (cornerx[i]-centerx)*ftemp;
	 cornery[i] = centery + (cornery[i]-centery)*ftemp;
	 }
   /* if any corner has x or y past limit, move the image */
   adjx = adjy = 0;
   for (i=0; i<4; ++i) {
      if (cornerx[i] > limit	 && (ftemp = cornerx[i] - limit) > adjx)
	 adjx = ftemp;
      if (cornerx[i] < 0.0-limit && (ftemp = cornerx[i] + limit) < adjx)
	 adjx = ftemp;
      if (cornery[i] > limit	 && (ftemp = cornery[i] - limit) > adjy)
	 adjy = ftemp;
      if (cornery[i] < 0.0-limit && (ftemp = cornery[i] + limit) < adjy)
	 adjy = ftemp;
      }
   if (calc_status == 2 && (adjx != 0 || adjy != 0) && (zwidth == 1.0))
      calc_status = 0;
   xxmin = cornerx[0] - adjx;
   xxmax = cornerx[1] - adjx;
   xx3rd = cornerx[2] - adjx;
   yymax = cornery[0] - adjy;
   yymin = cornery[1] - adjy;
   yy3rd = cornery[2] - adjy;
   adjust_corner(); /* make 3rd corner exact if very near other co-ords */
}

static void _fastcall smallest_add(double *num)
{
   *num += *num * 5.0e-16;
}

static int _fastcall ratio_bad(double actual, double desired)
{  double ftemp;
   if (desired != 0)
      if ((ftemp = actual / desired) < 0.95 || ftemp > 1.05)
	 return(1);
   return(0);
}


/* Save/resume stuff:

   Engines which aren't resumable can simply ignore all this.

   Before calling the (per_image,calctype) routines (engine), calcfract sets:
      "resuming" to 0 if new image, nonzero if resuming a partially done image
      "calc_status" to 1
   If an engine is interrupted and wants to be able to resume it must:
      store whatever status info it needs to be able to resume later
      set calc_status to 2 and return
   If subsequently called with resuming!=0, the engine must restore status
   info and continue from where it left off.

   Since the info required for resume can get rather large for some types,
   it is not stored directly in save_info.  Instead, memory is dynamically
   allocated as required, and stored in .fra files as a separate block.
   To save info for later resume, an engine routine can use:
      alloc_resume(maxsize,version)
	 Maxsize must be >= max bytes subsequently saved + 2; over-allocation
	 is harmless except for possibility of insufficient mem available;
	 undersize is not checked and probably causes serious misbehaviour.
	 Version is an arbitrary number so that subsequent revisions of the
	 engine can be made backward compatible.
	 Alloc_resume sets calc_status to 2 (resumable) if it succeeds; to 3
	 if it cannot allocate memory (and issues warning to user).
      put_resume({bytes,&argument,} ... 0)
	 Can be called as often as required to store the info.
	 Arguments must not be far addresses.
	 Is not protected against calls which use more space than allocated.
   To reload info when resuming, use:
      start_resume()
	 initializes and returns version number
      get_resume({bytes,&argument,} ... 0)
	 inverse of store_resume
      end_resume()
	 optional, frees the memory area sooner than would happen otherwise

   Example, save info:
      alloc_resume(sizeof(parmarray)+100,2);
      put_resume(sizeof(int),&row, sizeof(int),&col,
		 sizeof(parmarray),parmarray, 0);
    restore info:
      vsn = start_resume();
      get_resume(sizeof(int),&row, sizeof(int),&col, 0);
      if (vsn >= 2)
	 get_resume(sizeof(parmarray),parmarray,0);
      end_resume();

   Engines which allocate a large far memory chunk of their own might
   directly set resume_info, resume_len, calc_status to avoid doubling
   transient memory needs by using these routines.

   StandardFractal, calcmand, solidguess, and bound_trace_main are a related
   set of engines for escape-time fractals.  They use a common worklist
   structure for save/resume.  Fractals using these must specify calcmand
   or StandardFractal as the engine in fractalspecificinfo.
   Other engines don't get btm nor ssg, don't get off-axis symmetry nor
   panning (the worklist stuff), and are on their own for save/resume.

   */

#ifndef XFRACT
int put_resume(int len, ...)
#else
int put_resume(va_alist)
va_dcl
#endif
{
   va_list arg_marker;	/* variable arg list */
   char *source_ptr;
#ifdef XFRACT
   int len;
#endif

   if (resume_info == NULL)
      return(-1);
#ifndef XFRACT
   va_start(arg_marker,len);
#else
   va_start(arg_marker);
   len = va_arg(arg_marker,int);
#endif
   while (len)
   {
      source_ptr = va_arg(arg_marker,char *);
      far_memcpy(resume_info+resume_len,source_ptr,len);
      resume_len += len;
      len = va_arg(arg_marker,int);
   }
   return(0);
}

int alloc_resume(int alloclen, int version)
{
   if (resume_info != NULL) /* free the prior area if there is one */
      farmemfree(resume_info);
   if ((resume_info = farmemalloc((long)alloclen))== NULL)
   {
      static char msg[] = {"\
Warning - insufficient free memory to save status.\n\
You will not be able to resume calculating this image."};
      stopmsg(0,msg);
      calc_status = 3;
      return(-1);
   }
   resume_len = 0;
   put_resume(sizeof(int),&version,0);
   calc_status = 2;
   return(0);
}

#ifndef XFRACT
int get_resume(int len, ...)
#else
int get_resume(va_alist)
va_dcl
#endif
{
   va_list arg_marker;	/* variable arg list */
   char *dest_ptr;
#ifdef XFRACT
   int len;
#endif

   if (resume_info == NULL)
      return(-1);
#ifndef XFRACT
   va_start(arg_marker,len);
#else
   va_start(arg_marker);
   len = va_arg(arg_marker,int);
#endif
   while (len)
   {
      dest_ptr = va_arg(arg_marker,char *);
      far_memcpy(dest_ptr,resume_info+resume_offset,len);
      resume_offset += len;
      len = va_arg(arg_marker,int);
   }
   return(0);
}

int start_resume()
{
   int version;
   if (resume_info == NULL)
      return(-1);
   resume_offset = 0;
   get_resume(sizeof(int),&version,0);
   return(version);
}

void end_resume()
{
   if (resume_info != NULL) /* free the prior area if there is one */
   {
      farmemfree(resume_info);
      resume_info = NULL;
   }
}


/* Showing orbit requires converting real co-ords to screen co-ords.
   Define:
       Xs == xxmax-xx3rd	       Ys == yy3rd-yymax
       W  == xdots-1		       D  == ydots-1
   We know that:
       realx == lx0[col] + lx1[row]
       realy == ly0[row] + ly1[col]
       lx0[col] == (col/width) * Xs + xxmin
       lx1[row] == row * delxx
       ly0[row] == (row/D) * Ys + yymax
       ly1[col] == col * (0-delyy)
  so:
       realx == (col/W) * Xs + xxmin + row * delxx
       realy == (row/D) * Ys + yymax + col * (0-delyy)
  and therefore:
       row == (realx-xxmin - (col/W)*Xs) / Xv	 (1)
       col == (realy-yymax - (row/D)*Ys) / Yv	 (2)
  substitute (2) into (1) and solve for row:
       row == ((realx-xxmin)*(0-delyy2)*W*D - (realy-yymax)*Xs*D)
		      / ((0-delyy2)*W*delxx2*D-Ys*Xs)