|
asimes
|
 |
« on: March 11, 2012, 02:17:29 AM » |
|
I made my first Buddhabrot renderer. I just read some articles that explained it and then tried to make one, so I didn't really refer to anything technical. Here's a couple samples so far: "Finite" number of starting points (only uses the pixels on screen as random (x, y) points:  "Infinite" number of starting points (uses the pixels on the screen plus a random range):  Both of the images above were made by first removing all the pixels on the screen that are part of the Mandelbrot Set. After that's done they go to a second state for drawing. This state picks randomly among the remaining pixels, converts their location to an appropriate range (between -2 and 2), and draws their orbits. The second image adds a small random range to the points being picked to so that it has a bigger pool of numbers between -2, and 2. How do people normally handle this? Here is a the code of the random number picking: int randomPoint = nonMSet[int(random(nonMSetCount))];
int rx = randomPoint%width;
int ry = randomPoint/width;
float x = map(rx, 0, width, xmin, xmax);
float y = map(ry, 0, height, ymin, ymax);
x += random(-0.0001, 0.0001); // This is only if I want the second described point picking
y += random(-0.0001, 0.0001); // This is only if I want the second described point picking Another thing is that the pixels are all originally black (they have a value of 0). Whenever an orbit is tracked, it increments the appropriate pixel location (slowly increasing the value). A pixel that reaches a value of 255 is at pure white. Is there a method of automating when the process should be done? If I just leave it alone it will continue getting brighter until it is completely overblown. |
|
|
|
|
Logged
|
|
|
|
|
asimes
|
|
« Reply #1 on: March 11, 2012, 04:07:10 AM » |
|
Also discovered another problem... my method of finding the Mandelbrot Set was by testing every pixel to see if it never broke the bailout. Then I used the points that did break the bailout randomly when I did the drawing. The problem is that this makes slightly asymmetrical drawings (although it draws really fast). How do people normally avoid the Mandelbrot Set?
|
|
|
|
|
Logged
|
|
|
|
|
asimes
|
|
« Reply #2 on: March 11, 2012, 08:10:45 AM » |
|
Solved the point selection problem:
- Before the drawing process, test every pixel to see if it is part of the Mandelbrot Set
- Store every pixel location that is not part of the Mandelbrot Set in an array (I'll refer to this as nonMSet[])
- During the drawing process, randomly pick an element from nonMSet[] (the element is a point on the screen)
- Convert this point to the appropriate range for point testing (between -2 and 2 on the 'x' and 'y' axii)
- Add a random range to this point (so that all the points between pixels are considered)
- Increment the corresponding pixels[] element AND increment an array called undo[]
- If the point is found to be part of the Mandelbrot, set subtract the corresponding elements of undo[] from the matching elements of pixels[], otherwise do nothing pixels[]
It is almost as fast as the original method and tests every point between -2 and 2 on the 'x' and 'y' axii. If anyone is interested I will post code. Still interested in limiting the drawing from being too white / overexposed if anyone has advice. I'll put up a live version in a little while.
EDIT: Is it normal to need a little bit of "luck" to get a nice image with this process? Occasionally I get huge orbits at a time that kind of take over the image on high iterations.
|
|
|
|
« Last Edit: March 11, 2012, 08:16:55 AM by asimes »
|
Logged
|
|
|
|
Gluecker
Forums Newbie
Posts: 9
|
|
« Reply #3 on: March 11, 2012, 02:11:05 PM » |
|
I recommend (at least as long as you are not planning to make zooms):
[pseudocode]
1.) (variable type might be double, or doubledouble, etc)
ZRe = (random()-0.5)*4.0; //=> -2<Re<+2
ZIm = random()*2.0; //=> 0<Im<+2
2.)
if(isNotInMandelbrotSet(ZRe,ZIm)) //just a usual MSet EscapeTime with return type = boolean
{
//do MSet iteration here, each iteration do:
paint(ZRe,ZIm);
paint(ZRe,-ZIm); //(since it is symmetric)
}
3.) go to 1.)
What paint(re,im) does:
a) with clipping ( =the "pure white" )
Take an image and increase the brightness of the pixel corresponding to (re, im)
or:
b) Take a 2D-array (e.g. name it "map") with the dimensions equal to the image size, a variable counting the number of points and do:
pointsPainted++;
map[(int)(re*scale-ReOffset)][(int)(im*scale-ImOffset)]++;
now you can paint the array as a picture: Image.setPixelBrightness(x,y,map[x,y]/pointsPainted*maxBrightness);
obviously you will never reach a brightness > maxBrightness
I prefer to output the image in a separate thread. also, you can run 1) to 3) with many threads parallel.
you could also paint colors, depending on the iteration count, eg paint pixels blue if Iteration_count > x etc.
with random() in the code, you will need "luck" to get a nice image, but you could also just let the above code run for some hours, this should give you some nice details.
of course you will need to try out the parameters to get some good results (e.g. maxBrightness is usually greater than the max value the image can take, so if 255 is white, then maxBrightness is sth. like 500-1000 or so).
|
|
|
|
|
Logged
|
|
|
|
|
asimes
|
|
« Reply #4 on: March 11, 2012, 06:34:56 PM » |
|
When I said it maps to the appropriate range (-2, 2) I actually meant to the current zoom boundaries. I just thought it would be confusing to throw that in with the instructions. It does this (NOTE: 'dx' and 'dy' are floats with width and height of a pixel mapped to the complex plane): int randomPoint = nonMSet[int(random(nonMSetCount))];
int rx = randomPoint%width;
int ry = randomPoint/width;
float x = map(rx, 0, width, xmin, xmax)+random(-dx/2, dx/2);j
float y = map(ry, 0, height, ymin, ymax)+random(-dy/2, dy/2); The paint "mirroring" is a good idea, that definitely solves my problem of asymmetry. I like how the images look when I mix different max iterations with one another where each is a separate color. However, this is a post processing technique (Photoshop basically), how would someone do this while the program is running? |
|
|
|
|
Logged
|
|
|
|
|
asimes
|
|
« Reply #5 on: March 11, 2012, 10:21:49 PM » |
|
After making three 1,000,000 Iteration images, assigning each a RGB color, and mixing:  |
|
|
|
|
Logged
|
|
|
|
|
|
Gluecker
Forums Newbie
Posts: 9
|
|
« Reply #7 on: March 12, 2012, 01:49:01 AM » |
|
I like how the images look when I mix different max iterations with one another where each is a separate color. However, this is a post processing technique (Photoshop basically), how would someone do this while the program is running?
you could extend the 2D-Array to a 3D where the 3rd field has a size of 3 for red, green and blue. then you extend your equivalent of the paint(...) method by one parameter, the current iteration. => [...] paint(re,im,iter) { if(iter < 1000) map[converted_re][converted_im][0]++; //the converted_... is the mapping from the real and imaginary values to int else if(iter < 10000) map[converted_re][converted_im][1]++; else map[converted_re][converted_im][2]++; } Or you could use a double array instead of the int map and take a gaussian function with some color adjustment parameters. then add a small value to all of the 3 color fields depending on the iteration count. ( map[ x][y][0] += f(iter, param_red); etc, and come up with a clever f(...) ) good luck  |
|
|
|
|
Logged
|
|
|
|
ker2x
Fractal Molossus
Posts: 795
|
|
« Reply #8 on: March 13, 2012, 12:18:02 PM » |
|
Hummm.... One of the many possibilities to render a buddhabrot : 0) Create an "exposure" array (using screen dimension) 1) Select a random point from -2.0 to +2.0 (well, 2 random number, one for Real and one for Imaginary) 2) Check if it is in the mandelbrot set 3) If it is, go back to step 1. if it is not in the mandelbrot set (that's what you want), go to step 4. 4) Recompute again (from the very same point that you now know that it's not in the mandelbrot set) BUT : 4a) For each step (starting at iteration = 0 up to your Maximum Iteration you choosed) : 4b) convert the complex coordinate to screen coordinate but instead of drawing the point on the screen you increment by 1 the value (at the same coordinate) of your "exposure" array. 4c) stop when your point have a distance from (0,0) greater than 2. 5) Go back to step 1 and loop from 1 to 5 a few thousands times. (you don't want to update your screen for only 1 orbit, drawing is very slow) Now you have a 2D array full of values defining the number of times a point "landed" on a "pixel". If you draw this directly (you could) the screen will be very dark at the begining (the "exposure" value will be low) and too bright later (with value far above 255) So... 6) find what is the maximum value in your "exposure" array. then store this value (eg: in a variable called "maxExposure") 7) for each pixel on your screen check what is its exposure value (eg : in variable "pix") then convert it to a color. Assuming you're drawing with 255 grey level. that would be : pix/maxExposure*255  update the screen, go back to step 1 ( not step 0, you don't want to reset exposure values to 0). 9) Enjoy There is some problem using this colouring method. If for some reason a point received many many many exposure compared to the others, you will have a single very bright pixel (color 255) and everything else will be darkened (or completly black if the maxExposure is 255 times bigger than the 2nd brightest point). You may want to use some logaritmic formula (or other formula) to translate from "exposure" to color. Of course, it's all black and white for now. If you want colors : Create 3 exposure array instead of one : exposureR, exposureG, exposureB. You define a minIteration and maxIteration value for R, G, B : minIterR, maxIterR, minIterG, maxIterG, minIterB, maxIterB. eg : R : 1->100 G : 101 -> 200 B : 201 -> 300 You patch your step 4b and increment the one or array or the array depending of the value of the current iteration step. In this exemple, an orbit will start red, then become green, and will finish in blue. Not yet super nice, i highly suggest that each min/max Iter range overlap eachothers. Eg : R : 1 -> 100 G : 50 -> 150 B : 101 -> 200 (yes, at some point a single iteration may iterate exposure values of 2 or 3 arrays.) Now you can enjoy the result of this method : http://www.fractalforums.com/images-showcase-(rate-my-fractal)/ker2x's-buddhabrot-gallery/Hooooooo   |
|
|
|
|
Logged
|
|
|
|
ker2x
Fractal Molossus
Posts: 795
|
|
« Reply #9 on: March 13, 2012, 12:21:59 PM » |
|
You can do some easy early rejection to find if a point is in the mandelbrot set. My openCL code doing some test before trying the bruteforce method : bool isInMSet(
float cr,
float ci,
const uint maxIter,
const float escapeOrbit)
{
int iter = 0;
float zr = 0.0;
float zi = 0.0;
float ci2 = ci*ci;
float temp;
//Quick rejection check if c is in 2nd order period bulb
if( (cr+1.0) * (cr+1.0) + ci2 < 0.0625) return true;
//Quick rejection check if c is in main cardioid
float q = (cr-0.25)*(cr-0.25) + ci2;
if( q*(q+(cr-0.25)) < 0.25*ci2) return true;
// test for the smaller bulb left of the period-2 bulb
if (( ((cr+1.309)*(cr+1.309)) + ci*ci) < 0.00345) return true;
// check for the smaller bulbs on top and bottom of the cardioid
if ((((cr+0.125)*(cr+0.125)) + (ci-0.744)*(ci-0.744)) < 0.0088) return true;
if ((((cr+0.125)*(cr+0.125)) + (ci+0.744)*(ci+0.744)) < 0.0088) return true;
while( (iter < maxIter) && ((zr*zr+zi*zi) < escapeOrbit) )
{
temp = zr * zi;
zr = zr*zr - zi*zi + cr;
zi = temp + temp + ci;
iter++;
}
if( iter < maxIter)
{
return false;
} else {
return true;
}
} |
|
|
|
|
Logged
|
|
|
|
ker2x
Fractal Molossus
Posts: 795
|
|
« Reply #10 on: March 13, 2012, 12:27:34 PM » |
|
My first implementation of the method described above. (but it write in a .ppm file instead of drawing directly on screen. (the code between OPEN and CLOSE) In ... Fortran PROGRAM MANDEL
USE omp_lib
IMPLICIT NONE
CHARACTER ( len = 255 ), PARAMETER :: filename = 'buddhabrot.ppm'
INTEGER, PARAMETER :: file_out_unit = 10
INTEGER, PARAMETER :: n_max=100
INTEGER, PARAMETER :: grid_resolution = 512
INTEGER, PARAMETER :: grid_center = grid_resolution/2
INTEGER, PARAMETER :: zpower = 2
INTEGER, PARAMETER :: miniter = 2
INTEGER*8, PARAMETER :: batchSize = 10000000
REAL, PARAMETER :: escapeOrbit = 2
REAL, PARAMETER :: xmin = -1.0, xmax = 2.0, ymin = -1.3, ymax =1.3
REAL, PARAMETER :: intensityR = 2048.
REAL, PARAMETER :: intensityG = 2048.
REAL, PARAMETER :: intensityB = 2048.
!Track pixel exposure by color
INTEGER :: exposureRMap(grid_resolution, grid_resolution)
INTEGER :: exposureGMap(grid_resolution, grid_resolution)
INTEGER :: exposureBMap(grid_resolution, grid_resolution)
INTEGER :: maxRExposure, maxGExposure, maxBExposure
INTEGER :: minRExposure, minGExposure, minBExposure
INTEGER :: maxExposure, minExposure
INTEGER*8 :: i
REAL :: x,y
COMPLEX :: z, c
INTEGER :: iter
INTEGER :: tempX, tempY, tempYm
INTEGER :: ppm_i
INTEGER :: ppm_j
INTEGER :: ppm_jhi
INTEGER :: ppm_jlo
!initialize exposureMap to 1
exposureRMap = 0
exposureGMap = 0
exposureBMap = 0
!$OMP PARALLEL DO DEFAULT(PRIVATE),SHARED(exposureRMap, exposureGMap, exposureBMap)
DO i=1, batchSize
CALL RANDOM_NUMBER(x)
CALL RANDOM_NUMBER(y)
z = CMPLX(0,0)
!c = CMPLX(x*2.5 - 2 ,y*2.6 - 1.3) !choose a random point on complex plane
c = CMPLX((x*2.5 - 2) ,y*1.3 ) !choose a random point on complex plane
IF (notInMSet(c, n_max)) THEN !if it espace out of the mandelbrot set
DO iter=1, n_max !iterate and plot orbit
z = z**zpower + c !mandelbrot formula : Z = Z²+C
IF(iter .GE. miniter) THEN
TempX = INT(grid_resolution * (REAL(z) + xmax) / (xmax - xmin))
TempY = INT(grid_resolution * (AIMAG(z) + ymax) / (ymax - ymin))
TempYm = INT(grid_center - (TempY - grid_resolution/2))
IF((TempX > 0) .AND. (TempX < grid_resolution) .AND. (TempY > 0) .AND. (TempY < grid_resolution)) THEN
IF((iter > 2) .AND. (iter < 50)) THEN
exposureRMap(TempX, TempY) = exposureRMap(TempX, TempY) + 1
exposureRMap(TempX, TempYm) = exposureRMap(TempX, TempYm) + 1
END IF
IF((iter > 25) .AND. (iter < 75)) THEN
exposureGMap(TempX, TempY) = exposureGMap(TempX, TempY) + 1
exposureGMap(TempX, TempYm) = exposureGMap(TempX, TempYm) + 1
END IF
IF((iter > 50) .AND. (iter < 100)) THEN
exposureBMap(TempX, TempY) = exposureBMap(TempX, TempY) + 1
exposureBMap(TempX, TempYm) = exposureBMap(TempX, TempYm) + 1
ENDIF
END IF
END IF !(cabs(z)<4)
END DO
END IF
END DO
!$END PARALLEL
maxRExposure = MAXVAL(exposureRMap)
!minRExposure = MINVAL(exposureRMap)
maxGExposure = MAXVAL(exposureGMap)
!minGExposure = MINVAL(exposureGMap)
maxBExposure = MAXVAL(exposureBMap)
!minBExposure = MINVAL(exposureBMap)
!write(*,*) maxRExposure, minRExposure, maxGExposure, minGExposure, maxBExposure, minBExposure
!minExposure = MIN(minRExposure, minGExposure, minBExposure)
maxExposure = MAX(maxRExposure, maxGExposure, maxBExposure)
!exposureRMap = exposureRMap - minExposure
!exposureGMap = exposureGMap - minExposure
!exposureBMap = exposureBMap - minExposure
exposureRMap = (exposureRMap / REAL(maxRExposure))*intensityR
exposureGMap = (exposureGMap / REAL(maxGExposure))*intensityG
exposureBMap = (exposureBMap / REAL(maxBExposure))*intensityB
open ( unit = file_out_unit, file = filename, status = 'replace', &
form = 'formatted', access = 'sequential')
write ( file_out_unit, '(a2)' ) 'P3'
write ( file_out_unit, '(i5,2x,i5)' ) grid_resolution, grid_resolution
write ( file_out_unit, '(i5)' ) INT(MAX(intensityR,intensityG,intensityB))
do ppm_i = 1, grid_resolution
do ppm_jlo = 1, grid_resolution, 4
ppm_jhi = min ( ppm_jlo + 3, grid_resolution )
write ( file_out_unit, '(12i5)' ) &
( exposureRMap(ppm_i,ppm_j), exposureGMap(ppm_i,ppm_j),exposureBMap(ppm_i,ppm_j), ppm_j = ppm_jlo,ppm_jhi )
end do
end do
close ( unit = file_out_unit )
CONTAINS
PURE FUNCTION notInMset(c, n_max)
COMPLEX, INTENT(IN) :: c
INTEGER, INTENT(IN) :: n_max
INTEGER :: n
COMPLEX :: z
LOGICAL :: notInMSet
z = CMPLX(0,0)
n = 0
IF(((ABS(c - CMPLX(-1,0) )) < 0.25) .OR. ((ABS( 1.0 - SQRT(1-(4*c)) )) < 1.0 ) ) THEN
notInMset = .FALSE.
ELSE
DO WHILE (ABS(z) < escapeOrbit .AND. (n < n_max))
z = z**zpower + c
n = n + 1
END DO
IF (n >= n_max) THEN
notInMset = .FALSE.
ELSE
notInMset = .TRUE.
END IF
END IF
END FUNCTION notInMset
END |
|
|
|
|
Logged
|
|
|
|
|
asimes
|
|
« Reply #11 on: March 13, 2012, 04:31:12 PM » |
|
This process I posted above only finds about 0-6 points in the Mandelbrot Set on average for every 1,000 randomly points. It also doesn't have to double compute the orbit, it just subtracts from the pixels[] array if it turns out to be part of the Mandelbrot Set. - Before the drawing process, test every pixel to see if it is part of the Mandelbrot Set
- Store every pixel location that is not part of the Mandelbrot Set in an array (I'll refer to this as nonMSet[])
- During the drawing process, randomly pick an element from nonMSet[] (the element is a point on the screen)
- Convert this point to the appropriate range for point testing (between -2 and 2 on the 'x' and 'y' axii)
- Add a random range to this point (so that all the points between pixels are considered)
- Increment the corresponding pixels[] element AND increment an array called undo[]
- If the point is found to be part of the Mandelbrot, set subtract the corresponding elements of undo[] from the matching elements of pixels[], otherwise do nothing pixels[] I didn't have too much luck with the 3 array (for RGB) method unfortunately. Every time I got a big orbit drawn all at once it would be either exactly red, green, or blue and influenced the image too much for my liking. The big orbit sections tended to get overblown if I waited for them to also be drawn again by another color. I might try out the overlapping RGB array. If I set my max iterations to 1,000,000 then what would a good range be? The problem I'm seeing is that a lot of the orbits might not get anywhere near 1,000,000. If this is not a problem I was imagining: - Red range 0 to 500,000 - Green range 250,000 to 750,000 - Blue range 500,000 to 1,000,000 |
|
|
|
|
Logged
|
|
|
|
ker2x
Fractal Molossus
Posts: 795
|
|
« Reply #12 on: March 13, 2012, 05:00:36 PM » |
|
You're not looking for point that escape at exactly 1 Millions iterations but before 1 millions. You'll have a massive majority of orbit that will escape in the lower range. So with your setting above your picture will be 99% red You need to experiment, you may try something like this : 0->100.000 50.000-> 300.000 100.000 -> 1.000.000 No problem using very high iteration. I tried much higher using buddha++ (a buddhabrot renderer in Qt written by emilio and me, here : https://github.com/emmmile/buddha ) |
|
|
|
|
Logged
|
|
|
|
|
Non-parametric color mapping techniques
Programming
« 1 2 »
May 16, 2017, 12:55:21 AM
Dilber MC Theme by HarzeM