Mapowanie Fraktalne
W tym artykule postaram się przedstawić krok po kroku, jak za pomocą biblioteki Allegro, stworzyć prostą aplikację umożliwiającą rysowanie fraktali Mandelbrota, oraz jak nałożyć na fraktal własną bitmapę, tak by poddała się ona „fraktalnym zniekształceniom” :).
W omawianym programie nie kładłem specjalnego nacisku na jego skalowalność, ani przenośność. Skupiłem się tutaj jedynie na samej implementacji poszczególnych funkcji, nie na „ładnej” strukturze programu. Omawiany przykład będzie tworzony w środowisku Dev-C++, z wykorzystaniem wyżej wymienionej biblioteki Allegro. Opis konfiguracji środowiska do tworzenia aplikacji działających w trybie graficznym znajduje się w poprzednim artykule „Wstęp do Allegro”. Zatem do dzieła.
Tworzymy standardowy projekt Multimedia -> Allegro (static). Będziemy rozbudowywać nowo utworzony szablon. Dodajemy bibliotekę complex do programu (#include <complex>), która będzie odpowiedzialna za przetwarzanie liczb zespolonych.
W pierwszej kolejności musimy rozwiązać problem translacji współrzędnych. Wynika on z tego, że fraktal jest wykresem zbioru punktów na płaszczyźnie kartezjańskiej, w której początek układu jest w punkcie 0,0 (który dla nas jest jakby w środku pola widzenia), współrzędne każdego punktu są liczbami rzeczywistymi, wartości na osi x rosną od lewej do prawej, a wartości na osi y od dołu do góry. Ekran komputera natomiast, jest fragmentem płaszczyzny, której początek 0,0 jest umiejscowiony w górnym lewym rogu ekranu, współrzędne są liczbami całkowitymi z zakresów podyktowanych rozdzielczością aktualnego trybu graficznego. Dodatkowo oś y jest odwrócona, wartości w przeciwieństwie do układu kartezjańskiego rosną od góry do dołu. Potrzebne nam będą zatem dwie funkcje dokonujące translacji współrzędnych ekranu na współrzędne kartezjańskie. Zanim je napiszemy stworzymy jeszcze bardziej ogólną funkcje skalowania.
Funkcja skalowania ma za zadanie rzutować jedną wartość w pewnym układzie odniesienia na inną wartość w innym układzie odniesienia. Przykład:Wartość 10 w zakresie 0..100 rzutujemy na zakres 0..200. Nową wartością będzie 20. To jest coś w stylu przeliczania Celsjuszy na Fahrenheity :)
Aby napisać taką funkcje musimy znaleźć odpowiedni wzór. Mamy dane 2 pary liczb: dwie liczby opisujące pierwszy źródłowy przedział (A1,B1) i dwie liczby opisujące docelowy przedział (A2, B2 ), oraz jedną wartość X będącą pewną wartością w przedziale pierwszym. Musimy znaleźć taką funkcję y=aX+b, która dokona rzutowania X->Y z pierwszego przedziału na drugi przedział, musimy więc znaleźć współczynniki a i b. Rozwiązujemy następujący układ równań:
Otrzymujemy następujący wynik:
a = (B2-A2)/(B1-A1);
b=A2-aA1;
A więc funkcja skalująca będzie wyglądać tak:
inline double skaluj(double A1,double B1,double A2,double B2,double X)
{
double a=0,b=0;
if(B1==A1)return 0;
a=(B2-A2)/(B1-A1);
b=A2-a*A1;
return a*X+b;
}
Przed funkcją jest modyfikator inline ponieważ funkcja skaluj będzie wywoływana wiele tysięcy razy i nie zawiera długiego ciągu instrukcji. Dzięki temu funkcja zostanie „wklejona” w miejsca jej wywołań, zamiast „skakania” do niej. W ten sposób program się nieco wydłuża, ale za to zyskujemy na czasie, o czym będzie można się przekonać w trakcie generacji fraktala.
Teraz możemy zająć się funkcjami odpowiedzialnymi za translację współrzędnych ekranu na współrzędne kartezjańskie. Zanim jednak je napiszemy dokonamy zmian w dwóch miejscach programu. Za dyrektywą include umieszczamy dwie makrodefinicje:
#include <allegro.h> #include <complex> #define XRES 640 #define YRES 480
oraz zmieniamy argumenty w wywołaniu funkcji set_gfx_mode, z liczbowych na nasze zdefiniowane makrem:
res = set_gfx_mode(GFX_AUTODETECT_WINDOWED, XRES, YRES, 0, 0);
W ten sposób z jednego miejsca programu będziemy mogli zmienić rozdzielczość, bez konieczności ręcznego jej poprawiania w każdej funkcji odwołującej się do rozdzielczości. Stworzymy sobie teraz dwie funkcje translacji. Jedną z x ekranu na x kartezjańskie, i drugą z y ekranu na y kartezjańskie. W obu posłużymy się wcześniejszą funkcją skaluj.
inline double xs2xk(double x)
{
return skaluj(0,XRES,-1.0/zoom,1.0/zoom,x);
}
inline double ys2yk(double y)
{
return -skaluj(0,YRES,(-1.0/zoom)*ratio,(1.0/zoom)*ratio,y);
}
Brakuje dwóch obiektów: zoom i ratio. Zoom to powiększenie, a ratio to stosunek rozdzielczości pionowej do poziomej. Umieszczamy je w programie w przestrzeni globalnej:
const double ratio = double(YRES)/double(XRES); double zoom = 10.0;
Piszemy teraz funkcje testującą czy wszystko działa jak należy:
void test()
{
for(int la=0;la<XRES;la++)
for(int lb=0;lb<YRES;lb++)
{
double bright;
bright = sin(sqrt(pow(xs2xk(la),2)+pow(ys2yk(lb),2))*100.0);
bright = skaluj(-1.0,1.0,0,255,bright);
putpixel(screen,la,lb,makecol(bright,bright,bright));
}
}
I wywołujemy ją z tego miejsca programu w funkcji głownej:
init();
test();
while (!key[KEY_ESC]) {
/* put your code here */
}
Po uruchomieniu otrzymujemy poniższy obraz testowy:
Sprawdzamy czy zoom działa. Zmieniamy zoom z 10.0 na 2.0, wszystko działa:
Możemy teraz przystąpić do napisania zestawu funkcji odpowiedzialnych za generowanie naszego zbioru Mandelbrota. Popatrzmy na wzór matematyczny tego fraktala:
Jest to wzór rekurencyjny, który mówi, że każda kolejna liczba Z, jest równa poprzedniej liczbie Z podniesionej do kwadratu, zsumowanej z pewną liczbą C. Obie te liczby są liczbami zespolonymi. Zwykle początkowe Z jest równe 0, natomiast C zależy od współrzędnych badanego punktu. Każda liczba zespolona jest postaci a+ib, gdzie i jest jednostką urojoną, taką że i[sup]2[/sup]
= -1. Częścią rzeczywistą (Real) jest a, a częścią urojoną (Imaginary) b. Każdą liczbę zespoloną możemy przedstawić w postaci graficznej:
Odległość punktu Z od początku układu współrzędnych jest nazywana modułem liczby Z i oznaczana jako |Z|. Kąt w radianach pomiędzy osią rzeczywistą a półprostą przechodzącą przez punkt Z, mającą początek w początku układu współrzędnych nazywany jest argumentem liczby Z. Te dwa parametry będą miały dla nas duże znaczenie podczas nakładania tekstury na fraktal.
Wróćmy do naszego wzoru. Dla każdego piksela ekranu musimy wykonać pewną liczbę iteracji tego przypisania. Im więcej iteracji, tym bardziej szczegółowy będzie nasz obrazek, tym bardziej będziemy mogli się zagłębiać w nasz fraktal. Za Z podstawiamy na początku 0. Liczbę C tworzymy tak: jej część rzeczywista to współrzędna x, a jej część urojona to współrzędna y. Po kolejnych iteracjach moduł liczby może pozostać blisko zera, a może się zwiększyć bardzo gwałtowanie, stając się ogromną liczbą. Nas interesują przypadki w których moduł zostanie w przedziale od 0 do 2. Tylko gdy moduł będzie zawierał się w tym przedziale, będziemy mogli postawić na ekranie piksel.
Stwórzmy więc funkcję operującą na liczbach zespolonych, która dla podanych dwóch współrzędnych kartezjańskich zwróci finalną liczbę zespoloną Z. Zmienne zespolone tworzymy pisząc complex<precyzja>
nazwa_zmiennej; np. complex<int> liczba;
std::complex<double> oblicz(double x,double y)
{
std::complex<double> Z(0,0),C(x,y);
for(int la=0;la<iter;la++)
{
Z=Z*Z+C;
if(abs(Z)>cutoff)break;
}
return Z;
}
Odwołujemy się tutaj do liczby iteracji, oraz granicy odcięcia (cutoff). Stwórzmy więc zmienną opisującą liczbę iteracji oraz granicę odcięcia:
const double ratio = double(YRES)/double(XRES); double zoom = 0.8; int iter = 100; double cutoff = 2.0;
Zmieniliśmy również zoom, ponieważ fraktal będzie lepiej widoczny. Teraz możemy napisać główną funkcję rysującą fraktal, na razie będzie rysowała tylko białe lub czarne piksele:
void fractal()
{
clear(screen);
std::complex<double> wynik;
for(int la=0;la<XRES;la++)
for(int lb=0;lb<YRES;lb++)
{
wynik=oblicz(xs2xk(la),ys2yk(lb));
if(abs(wynik)<2.0)putpixel(screen,la,lb,makecol(255,255,255));
}
}
Na ekranie ukazuje się nam ciekawy kształt (po podstawieniu w miejsce wywołania funkcji test() funkcji fractal() ):
Ponieważ chcemy mieć możliwość wybierania sobie fragmentów fraktala do powiększenia, musimy stworzyć zmienne opisujące przesunięcie obrazu, oraz dokonać zmian w funkcjach odpowiedzialnych za translację współrzędnych.
const double ratio = double(YRES)/double(XRES); double zoom = 0.8; int iter = 100; double cutoff = 2.0; double xshift = 0.35; double yshift = 0.5;
Zmiany w kodzie funkcji translacji:
inline double xs2xk(double x)
{
return skaluj(0,XRES,-1.0/zoom,1.0/zoom,x)+xshift;
}
inline double ys2yk(double y)
{
return -skaluj(0,YRES,(-1.0/zoom)*ratio,(1.0/zoom)*ratio,y)-yshift;
}
Teraz możemy wybrać sobie dowolny fragment fraktala i zbliżyć się do niego:
Spróbujmy nieco ożywić nasz wykres :) Dodamy zielone, sinusoidalne przejście tonalne, wykorzystując drugi argument funkcji makecol (makecol(RED,GREEN,BLUE)), oraz skalowanie przedziału:
void fractal()
{
int bright;
double modul;
clear(screen);
std::complex<double> wynik;
for(int la=0;la<XRES;la++)
for(int lb=0;lb<YRES;lb++)
{
wynik=oblicz(xs2xk(la),ys2yk(lb));
modul=abs(wynik);
if(modul>cutoff)continue;
modul=sin(wynik.real()*6.0);
bright=skaluj(-1.0,1.0,0,255,modul);
putpixel(screen,la,lb,makecol(0,bright,0));
}
}
Generalnie, jeżeli chcemy uzyskać jakiś efekt graficzny, możemy uzależnić jedną z trzech składowych koloru od pewnego parametru, charakteryzującego dany punkt fraktala. Każdy punkt otrzymany z funkcji „oblicz” jest liczbą zespoloną. Możemy wydobyć z niej 4 różne parametry bezpośrednio, lub wykonać serię działań na tej liczbie i dopiero wydobyć parametry. Wykonanie serii działań matematycznych na tej liczbie zespolonej, może wprowadzić ciekawe efekty, które wpłyną na rozmieszczenie tekstury na naszym fraktalu. Cztery parametry, które możemy wydobyć z tej liczby to:
- część rzeczywista, zwracana przez wywołanie liczba.real();
- część urojona, zwracana przez wywołanie liczba.imag();
- jej moduł, otrzymywany przez przekazanie jej do funkcji abs(liczba);
- jej argument, otrzymywany przez przekazanie jej do funkcji arg(liczba);
Teraz spróbujemy okleić nasz twór teksturą. Jako tekstury użyjemy obrazu Monalisy:
Aby „okleić” powyższym obrazkiem nasz fraktal, musimy każdemu punktowi fraktala przyporządkować jakiś piksel z bitmapy. Ponieważ do określenia pozycji piksela potrzebne są dwie wartości, możemy określać pozycję piksela, który mamy zamiar wybrać poprzez odpowiednie dwa parametry punktu (liczby zespolonej) fraktala, np. moduł (amplitudę) i argument (fazę). Jedne z ładniejszych efektów uzyskuje się właśnie podczas przyporządkowywania amplitudowo-fazowego.
W pierwszej kolejności musimy zadbać o to, żeby nasza bitmapa została, załadowana. Proponuję użyć własnej bitmapy (np. własnego zdjęcia), oczywiście im wyższa rozdzielczość, tym lepiej. Kopiujemy plik bitmapy do katalogu projektu, nasz plik nazywa się monalisa.bmp.
Dodajemy do programu zmienną globalną typu BITMAP*, która będzie wskazywać na obszar pamięci przechowujący naszą teksturę (w praktyce jednak należy unikać stosowania obiektów globalnych, później są z nimi kłopoty, kod staje się mało przenośny, i przy rozbudowanym programie również czasochłonny pod względem analizy):
const double ratio = double(YRES)/double(XRES); double zoom = 10.8; int iter = 100; double cutoff = 2.0; double xshift = 0.32; double yshift = 0.5; BITMAP* tekstura;
Następnie w funkcji init() dodajemy wywołanie funkcji ładującej naszą bitmapę do pamięci:
install_mouse();
/* add other initializations here */
tekstura = load_bitmap("monalisa.bmp",0);
Teraz możemy napisać funkcję która zwraca wartość koloru piksela o podanych współrzędnych. Oczywiście można wywoływać zamiast tego funkcję getpixel, ale rozwiązanie to nie daje tylu możliwości co napisanie własnej funkcji ponieważ: - współrzędne możemy podawać za pomocą liczb zmiennoprzecinkowych, co da możliwość na późniejszą implementację wygładzania krawędzie pikseli za pomocą interpolacji (zostanie to omówione niebawem, posłużymy się interpolacją dwuliniową (bilinear) ) - takie „rozwarstwienie” programu ułatwi wprowadzenie efektów graficznych bazujących na barwach składowych koloru danego piksela, i wprowadzone zmiany w jednym miejscu będą powodowały skutki w całym oteksturowaniu fraktala. - możemy po swojemu określić zachowanie się tej funkcji w przypadku, kiedy współrzędne wykraczają poza obszar bitmapy. Np. czy ma wtedy zwracać czarny kolor, czy może kolor najbliższego możliwego piksela, czy przejść obrazek od drugiej strony, czy też może zastosować odbicia współrzędnych od krawędzi (czyli w momencie wyjścia poza obszar bitmapy, dalsze wychodzenie na zewnątrz, skutkuje cofaniem się)
A więc ta funkcja będzie wyglądać początkowo w ten sposób:
int getcolor(double x,double y)
{
int ix,iy;
//transformacja współrzędnych
ix = skaluj(-cutoff,+cutoff,0,tekstura->w,x);
iy = skaluj(-cutoff,+cutoff,0,tekstura->h,y);
//zwrocenie koloru
return getpixel(tekstura,ix,iy);
}
Spróbujemy jej użyć w funkcji fraktal, podstawiając za x część rzeczywistą, oraz za y część urojoną:
void fractal()
{
int kolor;
double a, b;
std::complex<double> wynik;
clear(screen);
for(int la=0;la<XRES;la++)
for(int lb=0;lb<YRES;lb++)
{
wynik=oblicz(xs2xk(la),ys2yk(lb));
a=wynik.real();
b=wynik.imag();
if(abs(wynik)>cutoff)continue;
kolor=getcolor(a,b);
putpixel(screen,la,lb,kolor);
}
}
Po zmniejszeniu liczby iteracji do 4, otrzymujemy taki obraz:
Jak widzimy metoda ta wygląda ciekawie. Rzuca się w oczy „pikselizacja” obrazu, której niebawem położymy kres. Natomiast jeszcze przez moment zajmijmy się funkcjami, odpowiedzialnymi za mapowanie tekstury. Przyjrzymy się teraz kilku sposobom mapowania.
Mapowanie nr 1
Metoda ta polega na „zablokowaniu” współrzędnej, gdy „próbuje” ona wyjść poza granice bitmapy. Kiedy współrzędna żądana będzie wynosić 300px, podczas gdy maksymalny wymiar obrazka wynosi 200px, współrzędna zwracana będzie równa 200px. Stąd smugi odchodzące od krawędzi obrazka.
Mapowanie nr 2
W tej metodzie mapowania, kiedy współrzędna próbuje wyjść za obszar bitmapy, skacze do przeciwległej krawędzi i jeżeli była rosnąca, to będzie rosła nadal, a jeżeli była malejąca, to będzie nadal maleć.
Mapowanie nr 3
To swoisty „ping-pong” :) Współrzędna, która próbuje przekroczyć granice bitmapy, zostaje odbita od krawędzi. Daje to złudzenie zniknięcia granic obrazków, nie widzimy ostrego odcięcia się jednych „płytek” obrazka od innych, tak jak to miało miejsce w poprzednim mapowaniu.
Mapowanie nr 4
Mapowanie to jest podobne do poprzedniego, ponieważ również współrzędna przekraczająca granicę obrazka zostaje odbita. Tutaj jednak występuje zależność nieliniowa między współrzędną żądaną, a zwracaną. Skutkiem tego jest deformacja obrazu, która nie występuje w poprzednich metodach. Zachowana jest tutaj ciągłość pierwszej pochodnej, a więc odbicia współrzędnych nie są „gwałtowne”. Obrazek obok jest stworzony na bazie sinusoidy.
W dalszej części będziemy używać wyłącznie mapowania „ping-pong”. Ogólnie rzecz biorąc, można tworzyć własne funkcje mapujące. Funkcje mapujące nie muszą być statyczne. Użycie dynamicznych funkcji mapujących opartych na szeregach Fouriera, może spowodować ciekawe efekty, na przykład podobne do falowania wody. Oto dwie funkcje które potrafią zamienić współrzędną żądaną na tą odpowiednią, według podanego numery trybu mapowania (mode: 0..3):
double xproc(double x,int mode)
{
switch(mode)
{
case 0 :
if(x>tekstura->w)x=tekstura->w-1;
if(x<0.0)x=0;
break;
case 1 :
if(x>0)x=fmod(x,tekstura->w);
else
{
x=fmod(-x,tekstura->w);
x=tekstura->w-1-x;
}
break;
case 2 :
if(x>0)x=fmod(x,tekstura->w*2);
else
{
x=fmod(-x,tekstura->w*2);
x=tekstura->w*2-2-x;
}
if(x>tekstura->w)x=-x+2*tekstura->w;
break;
case 3 :
x=sin(x/double(tekstura->w/2));
x=skaluj(-1.0,1.0,0,tekstura->w-1,x);
break;
}
return x;
}
I funkcja dla y:
double yproc(double y,int mode)
{
switch(mode)
{
case 0 :
if(y>tekstura->h)y=tekstura->h-1;
if(y<0.0)y=0;
break;
case 1 :
if(y>0)y=fmod(y,tekstura->h);
else
{
y=fmod(-y,tekstura->h);
y=tekstura->h-1-y;
}
break;
case 2 :
if(y>0)y=fmod(y,tekstura->h*2);
else
{
y=fmod(-y,tekstura->h*2);
y=tekstura->h*2-2-y;
}
if(y>tekstura->h)y=-y+2*tekstura->h;
break;
case 3 :
y=sin(y/double(tekstura->h/2));
y=skaluj(-1.0,1.0,0,tekstura->h-1,y);
break;
}
return y;
}
funkcja getcolor:
int getcolor(double x,double y)
{
double ix,iy;
//transformacja współrzędnych
ix = skaluj(-cutoff,+cutoff,0,tekstura->w,x);
iy = skaluj(-cutoff,+cutoff,0,tekstura->h,y);
ix=xproc(ix,2);
iy=yproc(iy,2);
//zwrocenie koloru
return getpixel(tekstura,int(ix),int(iy));
}
Teraz przyszedł czas na wygładzenie pikseli. Ponieważ współrzędna piksela jest na początku podawana zmiennoprzecinkowo, próby zbadania koloru piksela o współrzędnej rzeczywistej, będą zwracały kolor będący „mixem” kolorów sąsiednich pikseli o współrzędnych całkowitych. O tym do którego z sąsiednich kolorów będzie najbardziej podobny kolor piksela o którego pytamy, zadecyduje oczywiście odległość. Na poniższym schemacie przedstawiona jest metoda liniowej interpolacji kolorów, którą zastosujemy w projekcie:
Musimy napisać funkcję „interpoluj”, której argumentami będą współrzędne rzeczywiste i która będzie zwracać finalny kolor. Funkcja ta również będzie musiała sobie stworzyć cztery narożniki, których kolory będą znane. Posłużymy się tutaj matematycznymi funkcjami: podłoga (floor) i sufit (ceiling). Podłogą liczby 2.4 jest liczba 2.0 a sufitem liczba 3.0. Oto kod funkcji:
int interpoluj(double x,double y)
{
int k1,k2,k3,k4;//kolory
int r1,r2,r3,r4;//skladowe czerwone
int g1,g2,g3,g4;//skladowe zielone
int b1,b2,b3,b4;//skladowe niebieskie
int r,g,b,rl,gl,bl,rr,gr,br;//skladowe wynikowe, lewe i prawe
k2 = getpixel(tekstura,int(ceil(x)),int(floor(y)));
k1 = getpixel(tekstura,int(floor(x)),int(floor(y)));
k4 = getpixel(tekstura,int(floor(x)),int(ceil(y)));
k3 = getpixel(tekstura,int(ceil(x)),int(ceil(y)));
r1=getr(k1);
r2=getr(k2);
r3=getr(k3);
r4=getr(k4);
g1=getg(k1);
g2=getg(k2);
g3=getg(k3);
g4=getg(k4);
b1=getb(k1);
b2=getb(k2);
b3=getb(k3);
b4=getb(k4);
rl=skaluj(0,1.0,r2,r3,y-floor(y));
rr=skaluj(0,1.0,r1,r4,y-floor(y));
r=skaluj(0,1.0,rr,rl,x-floor(x));
gl=skaluj(0,1.0,g2,g3,y-floor(y));
gr=skaluj(0,1.0,g1,g4,y-floor(y));
g=skaluj(0,1.0,gr,gl,x-floor(x));
bl=skaluj(0,1.0,b2,b3,y-floor(y));
br=skaluj(0,1.0,b1,b4,y-floor(y));
b=skaluj(0,1.0,br,bl,x-floor(x));
return makecol(r,g,b);
}
Zmiany w funkcji getcolor:
int getcolor(double x,double y)
{
double ix,iy;
//transformacja współrzędnych
ix = skaluj(-cutoff,+cutoff,0,tekstura->w,x);
iy = skaluj(-cutoff,+cutoff,0,tekstura->h,y);
ix=xproc(ix,2);
iy=yproc(iy,2);
//zwrocenie koloru
return interpoluj(ix,iy);
}
Po wygenerowaniu fraktala, różnicę widać od razu:
I na tym koniec niniejszego artykułu. Zachęcam do własnoręcznego eksperymentowania z fraktalami i rozmaitymi efektami graficznymi. Czasami pisząc graficzny „freestyle”, możemy się natknąć na coś naprawdę ciekawego, nie mówiąc już o tym że jest to mile spędzony czas na wyrabianie w sobie odpowiedniego myślenia algorytmicznego. Dziękuję za uwagę :) Linki:
- http://komorra.110mb.com/fractaldream/ - Do pobrania program z GUI, oraz tapety.
- http://komorra.110mb.com/fractaldream/source.txt - Cały kod źródłowy omawianego programu.
- http://pl.wikipedia.org/wiki/Fraktal - Czym jest fraktal wg. Wikipedii.
- http://en.wikipedia.org/wiki/List_of_fractals_by_Hausdorff_dimension - Rodzaje fraktali.
- http://en.wikipedia.org/wiki/Beno%C3%AEt_Mandelbrot - Odkrywca opisywanej w tym artykule figury geometrycznej.
- http://mathworld.wolfram.com/ - Bardzo ciekawa encyklopedia matematyczna.
- http://komorra.110mb.com/ - Moja strona domowa.
- http://www.makzor.pl - Strona firmowa
Aby dodawać komentarze musisz być zalogowany!