/*************************************************************************** * * MandelVroom Generator * * Kevin L. Clague * * Copyright (C) 1987 * **************************************************************************/ #include "mand.h" #define PICTWINDXSIZE 100 #define PICTWINDYSIZE 100 #define PICTXSIZE (PICTWINDXSIZE - LEFTMARG - RIGHTMARG) #define PICTYSIZE (PICTWINDYSIZE - TOPMARG - BOTMARG) #define DEFAULTLEFT -2.0 #define DEFAULTRIGHT 2.0 #define DEFAULTTOP -2.0 #define DEFAULTBOT 2.0 extern struct Screen *screen; extern struct RastPort *rp; extern SHORT *CountBase; extern LONG NavTop, NavBot, NavLeft, NavRight; extern SHORT Zoom; extern struct Menu Menu, GenMenu; extern SHORT *ContourBase; extern UBYTE *ColorBase, MandType; struct Window *MandWind = (struct Window *) NULL; struct NewWindow NewMand = { 0,1, /* start position */ 80,80, /* width, height */ (UBYTE) 0, (UBYTE) 1, /* detail pen, block pen */ /* IDCMP flags */ MENUPICK | GADGETDOWN | GADGETUP | MOUSEBUTTONS, WINDOWDEPTH | WINDOWSIZING | ACTIVATE | REPORTMOUSE, /* MandWind flags */ (struct Gadget *) NULL, /* first gadget */ (struct Image *) NULL, /* user checkmark */ (UBYTE *) "Picture", /* MandWind title */ (struct Screen *) NULL, /* pointer to screen */ (struct BitMap *) NULL, /* pointer to superbitmap */ 80,80,320,200, /* sizing */ CUSTOMSCREEN /* type of screen */ }; float StartX = DEFAULTLEFT; float StartY = DEFAULTTOP; float EndX = DEFAULTRIGHT; float EndY = DEFAULTBOT; float GapY = ( DEFAULTBOT - DEFAULTTOP ) / PICTXSIZE; float GapX = 0.88 * ( DEFAULTBOT - DEFAULTTOP ) / PICTYSIZE; SHORT MaxCount = 1023; SHORT CountX = PICTXSIZE; SHORT CountY = PICTYSIZE; ULONG CalcTime; /* * Generate It */ GenerateIt() { ClosePalWind(); CloseContWind(); ClearMenuStrip(MandWind); SetMenuStrip(MandWind, &GenMenu); switch ( MandType ) { case 0: IntMandelbrot(); break; case 1: Mandelbrot(); break; } ClearMenuStrip(MandWind); SetMenuStrip(MandWind, &Menu); } /* * Fast Floating Point Mandelbrot Generator */ Mandelbrot() { int i, j, k, l; int TraceSize = 32; register float cura, curb, cura2, curb2; float curx, cury, napx, napy; float olda[64], oldb[64]; float *RealTrace,*ImagTrace; SHORT *CountPtr; UBYTE ColorXlate[1030]; #asm olda equ -548 oldb equ -292 TraceSize equ -20 l equ -16 k equ -12 curx equ -24 cury equ -28 #endasm /* Map contours into color translate table */ for (i = 0,j = 1029; i < NUMCONTS; i++) for (; j >= ContourBase[i] && j; ) ColorXlate[j--] = ColorBase[i]; while (j >= 0) ColorXlate[j--] = 0; /* free up old counts memory, get new picture size, get new counts memory */ GetCountsMemory(); /* clear the picture area */ SetAPen(MandWind->RPort, 0); RectFill(MandWind->RPort, LEFTMARG, TOPMARG, MandWind->Width - 12, MandWind->Height - 2); CountPtr = CountBase; /* calculate new picture's coordinates from zoom box */ ZoomIn(); /* start in the upper left hand corner */ cury = StartY; /* * for each pixel, calculate mandelbrot */ for (i = 0; i < CountY; i++) { curx = StartX; for (j = 0; j < CountX; j++) { cura = curb = cura2 = curb2 = 0; for (k = 0; k < MaxCount; k += TraceSize) { #ifdef SLOW RealTrace = olda; ImagTrace = oldb; for (l = 0; l < TraceSize; l++) { *(RealTrace++) = cura; *(ImagTrace++) = curb; curb *= cura; curb += curb + cury; cura = cura2 - curb2 + curx; cura2 = cura * cura; curb2 = curb * curb; if (cura2+curb2 > 4.0) goto out; } RealTrace = olda; ImagTrace = oldb; for (l = 0; l < TraceSize; l++) if (cura == *(RealTrace++) && curb == *(ImagTrace++)) { k += MaxCount; goto out; } #else #asm public _LVOSPMul public _LVOSPAdd public _LVOSPSub public _LVOSPCmp public _MathBase move.l _MathBase,a6 lea.l olda(a5),a0 ; RTrace = olda lea.l oldb(a5),a1 ; ITrace = oldb clr.l l(a5) ;for (l = TraceSize; l--; ) { lloop move.l d4,(a0)+ ; *(RTrace++) = cura move.l d5,(a1)+ ; *(ITrace++) = curb move.l d4,d0 ; curb *= cura move.l d5,d1 jsr _LVOSPMul(a6) move.l d0,d1 ; curb += curb jsr _LVOSPAdd(a6) move.l cury(a5),d1 ; curb += cury jsr _LVOSPAdd(a6) move.l d0,d5 move.l d6,d0 ; cura = cura2-curb2+curx move.l d7,d1 jsr _LVOSPSub(a6) move.l curx(a5),d1 jsr _LVOSPAdd(a6) move.l d0,d4 ; move.l d0,d1 ; cura2 = cura*cura jsr _LVOSPMul(a6) move.l d0,d6 move.l d5,d0 ; curb2 = curb*curb move.l d5,d1 jsr _LVOSPMul(a6) move.l d0,d7 move.l d6,d1 ; if (cura2+curb2 > 4.0) goto out: jsr _LVOSPAdd(a6) move.l #$80000043,d1 jsr _LVOSPCmp(a6) bgt out add.l #1,l(a5) move.l l(a5),d0 cmp.l TraceSize(a5),d0 blt lloop ; lea olda(a5),a0 ; RTrace = olda lea oldb(a5),a1 ; ITrace = oldb clr.l l(a5) ; for (l = TraceSize; l--; ) lloop1 cmp.l (a0)+,d4 bne skipit cmp.l (a1)+,d5 bne nextl move.w _MaxCount,d0 ext.l d0 move.l d0,k(a5) move.l #0,l(a5) bra out skipit add.l #4,a1 nextl add.l #1,l(a5) move.l l(a5),d0 cmp.l TraceSize(a5),d0 blt lloop1 #endasm #endif } k = 0; l = MaxCount; out: #asm out #endasm SetAPen(rp, ColorXlate[ k + l ]); WritePixel(rp, j + LEFTMARG,i + TOPMARG); if (CountPtr) *(CountPtr++) = k + l; curx += GapX; if ( CheckForStop() ) { DisplayBeep(screen); return(0); } } cury += GapY; } DisplayBeep(screen); } /* Mandelbrot */ /* * 32 bit fixed point generator */ IntMandelbrot() { LONG i, j, k; LONG curx, cury, gapx, gapy; register LONG cura, curb, cura2, curb2; float gap; SHORT *CountPtr; LONG OldA[32],OldB[32]; UBYTE MandFlag; UBYTE ColorXlate[1030]; UBYTE OldPen, NewPen; LONG xl,yl; OldPen = 0xff; #define BITS2SHIFT 27 /* Map contours into color translate table */ for (i = 0,j = 1029; i < NUMCONTS; i++) for (; j >= ContourBase[i] && j; ) ColorXlate[j--] = ColorBase[i]; while (j >= 0) ColorXlate[j--] = 0; /* free up old counts memory, get new picture size, get new counts memory */ GetCountsMemory(); /* clear the picture area */ SetAPen(MandWind->RPort, 0); RectFill(MandWind->RPort, LEFTMARG, TOPMARG, MandWind->Width - 12, MandWind->Height - 2); CountPtr = CountBase; /* calculate new picture's coordinates from zoom box */ ZoomIn(); /* figure out horizontal and verticle distances between points in plane. * convert them to fixed point format */ gapy = (int) (GapY*((float)(1<> 12; move.l d1,d2 ; get copy of op1 sign bit bpl pos1 ; get absolute value of op1 neg.l d1 pos1 eor.l d5,d2 ; calculate result sign tst.l d5 ; get absolute value of op2 bpl pos2 neg.l d5 pos2 move.l d1,d0 ; get a copy of op1 swap d0 ; get high half of op1 move.w d0,d2 ; save a copy of high half mulu d5,d0 ; multiply op2 low by op1 high clr.w d0 ; clear least significant part swap d0 ; put it in it's place swap d5 ; get high half of op2 mulu d5,d1 ; multiply op2 high with op1 low clr.w d1 ; clear least significant part swap d1 ; put it in its place mulu d2,d5 ; multiply op2 high by op1 high add.l d0,d5 ; add partial results add.l d1,d5 ; add partial results tst.l d2 ; is the result negative? bpl pos3 neg.l d5 ; yes, better negate it. pos3 asl.l d6,d5 ; now, rescale it. ; ; curb += curb + cury; add.l d5,d5 ; double it and add cury add.l -20(a5),d5 ; ; cura2 = cura * cura; move.l d4,d0 ; get absolute value of a in d0 bpl posa neg.l d0 posa move.l d0,d1 ; copy absolute value into d1 swap d1 ; get high part in d1 mulu d1,d0 ; multiply high and low destroying low clr.w d0 ; clear the least significant part swap d0 ; put most sig. part in low half mulu d1,d1 ; multiply high and high destroing high add.l d0,d1 ; add in lower half twice add.l d0,d1 asl.l d6,d1 ; get radix point back in correct place bvs bailout ; ; curb2 = curb * curb; move.l d5,d0 ; get absolute value of a in d0 bpl posb neg.l d0 posb move.l d0,d2 ; copy absolute value into d2 swap d2 ; get high part in d2 mulu d2,d0 ; multiply high and low destroying low clr.w d0 ; clear the least significant part swap d0 ; put most sig. part in low half mulu d2,d2 ; multiply high and high destroing high add.l d0,d2 ; add in lower half twice add.l d0,d2 asl.l d6,d2 ; get radix point back in correct place bvs bailout ; move.l d1,d0 ; if (cura2 + curb2 >= 4) goto bailout; add.l d2,d0 bvs bailout cmp.l #536870912,d0 bge bailout dbra d7,KLoop ; keep going ; move.w #Trace,d7 lea OldA(a5),a0 ; lea OldB(a5),a1 SLoop cmp.l (a0)+,d4 ; check A ; bne NotIt ; cmp.l (a1),d5 ; check B beq GotIt NotIt add #2,a1 dbra d7,SLoop ; check next one add.w #Trace,k(a5) ; drop the iteration count cmp.w #1023-Trace,k(a5) blt AddrLoop ; GotIt move.w #1023,d0 bra Adjusted bailout move.b #0,-293(a5) ; clear MandFlag here move.w #Trace,d0 sub.w d7,d0 add.w k(a5),d0 Adjusted ext.l d0 move.l d0,-12(a5) #endasm ;} else { /* This code does mandelbrot without any trace lookup. * It is the fastest generator in the house. */ #asm ; ; d1 - a2 ; d2 - b2 ; d4 - a ; d5 - b ; d6 - BITS2SHIFT ; d7 - k ; ; cura = curb = curc = curd = 0; moveq #0,d1 moveq #0,d2 moveq #0,d4 moveq #0,d5 move.l #Bits2Shift,d6 move.w #1023,d7 FKLoop ; ; cura = cura2 - curb2 + curx; exg d1,d4 ; exchange cura and cura2 sub.l d2,d4 ; subtract curb add.l -16(a5),d4 ; add curx ; ; curb = cura * curb >> 12; move.l d1,d2 ; get copy of op1 sign bit bpl Fpos1 ; get absolute value of op1 neg.l d1 Fpos1 eor.l d5,d2 ; calculate result sign tst.l d5 ; get absolute value of op2 bpl Fpos2 neg.l d5 Fpos2 move.l d1,d0 ; get a copy of op1 swap d0 ; get high half of op1 move.w d0,d2 ; save a copy of high half mulu d5,d0 ; multiply op2 low by op1 high clr.w d0 ; clear least significant part swap d0 ; put it in it's place swap d5 ; get high half of op2 mulu d5,d1 ; multiply op2 high with op1 low clr.w d1 ; clear least significant part swap d1 ; put it in its place mulu d2,d5 ; multiply op2 high by op1 high add.l d0,d5 ; add partial results add.l d1,d5 ; add partial results tst.l d2 ; is the result negative? bpl Fpos3 neg.l d5 ; yes, better negate it. Fpos3 asl.l d6,d5 ; now, rescale it. ; ; curb += curb + cury; add.l d5,d5 ; double it and add cury add.l -20(a5),d5 ; ; cura2 = cura * cura; move.l d4,d0 ; get absolute value of a in d0 bpl Fposa neg.l d0 Fposa move.l d0,d1 ; copy absolute value into d1 swap d1 ; get high part in d1 mulu d1,d0 ; multiply high and low destroying low clr.w d0 ; clear the least significant part swap d0 ; put most sig. part in low half mulu d1,d1 ; multiply high and high destroing high add.l d0,d1 ; add in lower half twice add.l d0,d1 asl.l d6,d1 ; get radix point back in correct place bvs Fbailout ; ; curb2 = curb * curb; move.l d5,d0 ; get absolute value of a in d0 bpl Fposb neg.l d0 Fposb move.l d0,d2 ; copy absolute value into d2 swap d2 ; get high part in d2 mulu d2,d0 ; multiply high and low destroying low clr.w d0 ; clear the least significant part swap d0 ; put most sig. part in low half mulu d2,d2 ; multiply high and high destroing high add.l d0,d2 ; add in lower half twice add.l d0,d2 asl.l d6,d2 ; get radix point back in correct place bvs Fbailout ; move.l d1,d0 ; if (cura2 + curb2 >= 4) goto bailout; add.l d2,d0 bvs Fbailout cmp.l #536870912,d0 bge Fbailout dbra d7,FKLoop move.b #1,-293(a5) ; set MandFlag here Fbailout move.w #1023,d0 sub.w d7,d0 move.w d0,d7 ext.l d7 move.l d7,-12(a5) #endasm } NewPen = ColorXlate[k]; if (NewPen != OldPen) { SetAPen(rp,NewPen); OldPen = NewPen; } if (NewPen) WritePixel(rp, j, i); if (CountPtr) *(CountPtr++) = k; curx += gapx; if ( CheckForStop() ) { /* Clean up here */ DisplayBeep(screen); return(0); } } cury += gapy; } DisplayBeep(screen); } /* Mandelbrot */ /* * Free up old memory and get memory for new picture */ GetCountsMemory() { /* free up old pictures iteration count pile */ if (CountBase) { FreeMem(CountBase,CountX*CountY*sizeof(SHORT)); CountBase = (SHORT *) NULL; } /* figure out new picture size */ CountX = MandWind->Width-(LEFTMARG+RIGHTMARG); CountY = MandWind->Height-(TOPMARG+BOTMARG); /* Allocate memory for new picture */ CountBase = (SHORT *) AllocMem(CountX*CountY*sizeof(SHORT),MEMF_CLEAR); if (CountBase == (SHORT *) NULL) { DispErrMsg("Can't save counts. Out of RAM!!",0); return(1); } return(0); } /* * ZoomIn */ ZoomIn() { if (Zoom) { EndY = StartY + GapY * (float)(NavBot-TOPMARG); EndX = StartX + GapX * (float)(NavRight-LEFTMARG); StartY = StartY + GapY * (float)(NavTop-TOPMARG); StartX = StartX + GapX * (float)(NavLeft-LEFTMARG); Zoom = 0; } CalculateGaps(); } /* ZoomIn */ /* * Calculate Gaps */ CalculateGaps() { GapY = (EndY-StartY)/(float)CountY; if (screen->ViewPort.Modes & HIRES) if (screen->ViewPort.Modes & INTERLACE) GapX = GapY*0.88; else GapX = GapY*0.44; else if (screen->ViewPort.Modes & INTERLACE) GapX = GapY*1.76; else GapX = GapY*0.88; } /* * Check For Stop */ CheckForStop() { struct IntuiMessage *message; struct Gadget *gadget; ULONG class; USHORT code; message = (struct IntuiMessage *) GetMsg(MandWind->UserPort); if (message) { class = message->Class; code = message->Code; ReplyMsg(message); if (class == MENUPICK && code != MENUNULL) { ClearMenuStrip(MandWind); return(1); } } return(0); } /* Check for stop */ /* * Open the Mandelbrot window */ OpenMandWind() { extern struct Window *OpenMyWind(); if (MandWind == (struct Window *) NULL) { MandWind = OpenMyWind(&NewMand, screen, NULL, 109, 78, 320, 200); SetMenuStrip(MandWind, &Menu); ForceNormPointer(); } } /* OpenMandWind */ /* * Close the Mandelbrot Window */ CloseMandWind() { if (MandWind != (struct Window *) NULL) { CloseMyWind(MandWind, NewMand.FirstGadget); MandWind = (struct Window *) NULL; } } /* CloseMandWind */