找不到按照我的目标为 Mandelbrot 上色的方法
Posted
技术标签:
【中文标题】找不到按照我的目标为 Mandelbrot 上色的方法【英文标题】：Can't find a way to color the Mandelbrot-set the way i'm aiming for 【发布时间】：2019-10-05 08:04:16 【问题描述】：我已经成功地为 Mandelbrot 集着色，尽管我无法放大很远，直到它变得“模糊”并且图案停止。我通过增加 max_iteration 来解决这个问题，这可行，但是在 *1 放大倍率时我得到的颜色很少，而且只有在放大时才会出现很多颜色。我理解为什么会发生这种情况，因为在“真正的”Mandelbrot-set 中没有颜色和增加 max_iterations 只会让它更接近这一点。但我的问题是，像 youtube 这样的缩放如何在整个缩放过程中拥有美丽的色彩，同时仍然能够放大到永远的感觉？
我尝试在网上到处寻找，但找不到解决方案，当我查看这些 youtube 缩放的描述时，他们似乎几乎没有提供有关他们如何进行缩放的任何信息。
这只是绘制 Mandelbrot 集的代码部分。下面的代码是在处理过程中编写的，它是添加了可视化库的 java。您可以在此处了解有关该计划的更多信息：https://processing.org/
//m is max_iterations
//mb is the function which calculates how many iterations each point took to escape to infinity. I won't be including the function since I know it works fine and it's quite messy.
//i'm using a HSB/HSV system to draw the mandelbrot
hue=(mb(x, y, m)*360)/m;
sat=255;
if (mb(x, y, m)<m)
val=255;
else
val=0;
stroke(hue,sat,val);
point(x, y);
我明白为什么会出现我的问题，但我不知道如何解决。
这是一个具有低 max_iterations 并缩小的图像，您可以看到它非常丰富多彩：
这是一个 max_iterations 较低的图像，稍微放大了一点，你可以看到它很无聊而且不是很丰富多彩：
这是一个具有高 max_iterations 并缩小的图像，您可以看到它不是很彩色：
这是一个具有高 max_iterations 并放大的图像，您可以看到它非常丰富多彩：
【问题讨论】：
你能发一张问题的图片吗？ 还创建并发布您的问题的minimal reproducible example，类似于我的Java Swing Mandelbrot MCVE。 这是我计算 Mandelbrot 集颜色的方法：Color.getHSBColor(i / 255F, 1, 1) 其中 i 是发散的迭代次数。我可能还建议在放大时动态增加迭代次数。
我现在正在添加图像。
为什么不限制色调变化，而是改变饱和度和亮度
【参考方案1】：
先看看这个相关的QA：
Mandelbrot Set - Color Spectrum Suggestions?主要思想是使用直方图将颜色渐变更有效地分配给已使用的索引，而不是在未使用的索引上统一浪费许多颜色。它还使用了一个特定的视觉上令人愉悦的渐变函数：
RGB values of visible spectrum 其他人建议的动态最大迭代次数只会影响整体性能和缩放细节。但是，如果您想要没有缩放的漂亮颜色，那么您需要计算 浮点迭代计数，这也称为 Mandelbrot Escape。有一种数学方法可以从方程的最后一个子结果计算迭代计数小数部分。欲了解更多信息，请参阅：
Renormalizing the Mandelbrot Escape但是我从来没有尝试过，所以带着偏见阅读：如果我没看错，你想要计算这个方程：
mu = m + frac = n + 1 - log (log |Z(n)|) / log 2
n 是您的迭代次数，Z(n) 是您正在迭代的方程的复杂域子结果。所以现在从 mu 计算颜色，现在是浮点数，而不是从 n...
[Edit2] GLSL mandelbrot 基于上面的链接进行分数转义
我添加了分数转义并修改了直方图多通道重新着色以匹配新输出...
顶点：
// Vertex
#version 420 core
layout(location=0) in vec2 pos; // glVertex2f <-1,+1>
out smooth vec2 p; // texture end point <0,1>
void main()
p=pos;
gl_Position=vec4(pos,0.0,1.0);
片段：
// Fragment
#version 420 core
uniform vec2 p0=vec2(0.0,0.0); // mouse position <-1,+1>
uniform float zoom=1.000; // zoom [-]
uniform int n=100; // iterations [-]
uniform int sh=7; // fixed point accuracy [bits]
uniform int multipass=0; // multi pass?
in smooth vec2 p;
out vec4 col;
const int n0=1; // forced iterations after escape to improve precision
vec3 spectral_color(float l) // RGB <0,1> <- lambda l <400,700> [nm]
float t; vec3 c=vec3(0.0,0.0,0.0);
if ((l>=400.0)&&(l<410.0)) t=(l-400.0)/(410.0-400.0); c.r= +(0.33*t)-(0.20*t*t);
else if ((l>=410.0)&&(l<475.0)) t=(l-410.0)/(475.0-410.0); c.r=0.14 -(0.13*t*t);
else if ((l>=545.0)&&(l<595.0)) t=(l-545.0)/(595.0-545.0); c.r= +(1.98*t)-( t*t);
else if ((l>=595.0)&&(l<650.0)) t=(l-595.0)/(650.0-595.0); c.r=0.98+(0.06*t)-(0.40*t*t);
else if ((l>=650.0)&&(l<700.0)) t=(l-650.0)/(700.0-650.0); c.r=0.65-(0.84*t)+(0.20*t*t);
if ((l>=415.0)&&(l<475.0)) t=(l-415.0)/(475.0-415.0); c.g= +(0.80*t*t);
else if ((l>=475.0)&&(l<590.0)) t=(l-475.0)/(590.0-475.0); c.g=0.8 +(0.76*t)-(0.80*t*t);
else if ((l>=585.0)&&(l<639.0)) t=(l-585.0)/(639.0-585.0); c.g=0.84-(0.84*t) ;
if ((l>=400.0)&&(l<475.0)) t=(l-400.0)/(475.0-400.0); c.b= +(2.20*t)-(1.50*t*t);
else if ((l>=475.0)&&(l<560.0)) t=(l-475.0)/(560.0-475.0); c.b=0.7 -( t)+(0.30*t*t);
return c;
void main()
int i,j,N;
vec2 pp;
float x,y,q,xx,yy,mu;
pp=(p/zoom)-p0; // y (-1.0, 1.0)
pp.x-=0.5; // x (-1.5, 0.5)
for (x=0.0,y=0.0,xx=0.0,yy=0.0,i=0;(i<n-n0)&&(xx+yy<4.0);i++)
q=xx-yy+pp.x;
y=(2.0*x*y)+pp.y;
x=q;
xx=x*x;
yy=y*y;
for (j=0;j<n0;j++,i++) // 2 more iterations to diminish fraction escape error
q=xx-yy+pp.x;
y=(2.0*x*y)+pp.y;
x=q;
xx=x*x;
yy=y*y;
mu=float(i)-log(log(sqrt(xx+yy))/log(2.0));
mu*=float(1<<sh); i=int(mu);
N=n<<sh;
if (i>N) i=N;
if (i<0) i=0;
if (multipass!=0)
// i
float r,g,b;
r= i &255; r/=255.0;
g=(i>> 8)&255; g/=255.0;
b=(i>>16)&255; b/=255.0;
col=vec4(r,g,b,255);
else
// RGB
q=float(i)/float(N);
q=pow(q,0.2);
col=vec4(spectral_color(400.0+(300.0*q)),1.0);
CPU 端 C++/VCL 代码：
//---------------------------------------------------------------------------
#include <vcl.h>
#pragma hdrstop
#include "Unit1.h"
#include "gl\\OpenGL3D_double.cpp"
//---------------------------------------------------------------------------
#pragma package(smart_init)
#pragma resource "*.dfm"
TForm1 *Form1;
OpenGLscreen scr;
GLSLprogram shd;
float mx=0.0,my=0.0,mx0=0.0,my0=0.0,mx1=0.0,my1=0.0;
TShiftState sh0,sh1;
int xs=1,ys=1;
float zoom=1.000;
int sh=7;
int N=256;
int _multi=0;
unsigned int queryID[2];
#define multi_pass
OpenGLtexture txr;
//---------------------------------------------------------------------------
DWORD spectral_color(float l) // RGB <0,1> <- lambda l <400,700> [nm]
float t; float r,g,b; DWORD c,x; r=0.0; g=0.0; b=0.0;
if ((l>=400.0)&&(l<410.0)) t=(l-400.0)/(410.0-400.0); r= +(0.33*t)-(0.20*t*t);
else if ((l>=410.0)&&(l<475.0)) t=(l-410.0)/(475.0-410.0); r=0.14 -(0.13*t*t);
else if ((l>=545.0)&&(l<595.0)) t=(l-545.0)/(595.0-545.0); r= +(1.98*t)-( t*t);
else if ((l>=595.0)&&(l<650.0)) t=(l-595.0)/(650.0-595.0); r=0.98+(0.06*t)-(0.40*t*t);
else if ((l>=650.0)&&(l<700.0)) t=(l-650.0)/(700.0-650.0); r=0.65-(0.84*t)+(0.20*t*t);
if ((l>=415.0)&&(l<475.0)) t=(l-415.0)/(475.0-415.0); g= +(0.80*t*t);
else if ((l>=475.0)&&(l<590.0)) t=(l-475.0)/(590.0-475.0); g=0.8 +(0.76*t)-(0.80*t*t);
else if ((l>=585.0)&&(l<639.0)) t=(l-585.0)/(639.0-585.0); g=0.84-(0.84*t) ;
if ((l>=400.0)&&(l<475.0)) t=(l-400.0)/(475.0-400.0); b= +(2.20*t)-(1.50*t*t);
else if ((l>=475.0)&&(l<560.0)) t=(l-475.0)/(560.0-475.0); b=0.7 -( t)+(0.30*t*t);
r*=255.0; g*=255.0; b*=255.0;
x=r; c =x;
x=g; c|=x<<8;
x=b; c|=x<<16;
return c;
//---------------------------------------------------------------------------
void gl_draw()
scr.cls();
// matrix for old GL rendering
glMatrixMode(GL_PROJECTION);
glLoadIdentity();
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
glMatrixMode(GL_TEXTURE);
glLoadIdentity();
// GLSL uniforms
shd.bind();
shd.set2f("p0",mx,my); // pan position
shd.set1f("zoom",zoom); // zoom
shd.set1i("n",N); // iterations
shd.set1i("sh",sh); // fixed point accuracy (shift)
shd.set1i("multipass",_multi); // single/multi pass
// issue the first query
// Records the time only after all previous
// commands have been completed
glQueryCounter(queryID[0], GL_TIMESTAMP);
// QUAD covering screen
glColor3f(1.0,1.0,1.0);
glBegin(GL_QUADS);
glVertex2f(-1.0,+1.0);
glVertex2f(-1.0,-1.0);
glVertex2f(+1.0,-1.0);
glVertex2f(+1.0,+1.0);
glEnd();
shd.unbind();
// [multipas]
if (_multi)
float t,m,n=N<<sh;
DWORD *hist=new DWORD[n+1];
int sz=txr.xs*txr.ys,i,j;
// get rendered image
glReadPixels(0,0,txr.xs,txr.ys,GL_RGBA,GL_UNSIGNED_BYTE,txr.txr);
// compute histogram
for (i=0;i<=n;i++) hist[i]=0;
for (i=0;i<sz;i++) hist[txr.txr[i]&0x00FFFFFF]++;
// histogram -> used color index (skip holes)
for (i=1,j=1;i<=n;i++)
if (hist[i]) hist[i]=j; j++;
// used color index -> color
m=1.0/float(j); hist[0]=0x00000000;
for (i=1;i<=n;i++)
if (hist[i]) t=hist[i]; t*=m; hist[i]=spectral_color(400.0+(300.0*t));
else hist[i]=0x00000000;
// recolor image
for (i=0;i<sz;i++) txr.txr[i]=hist[txr.txr[i]&0x00FFFFFF];
// render it back
scr.cls();
txr.bind();
glColor3f(1.0,1.0,1.0);
glBegin(GL_QUADS);
glTexCoord2f(0.0,1.0); glVertex2f(-1.0,+1.0);
glTexCoord2f(0.0,0.0); glVertex2f(-1.0,-1.0);
glTexCoord2f(1.0,0.0); glVertex2f(+1.0,-1.0);
glTexCoord2f(1.0,1.0); glVertex2f(+1.0,+1.0);
glEnd();
txr.unbind();
glDisable(GL_TEXTURE_2D);
delete[] hist;
// issue the second query
// records the time when the sequence of OpenGL
// commands has been fully executed
glQueryCounter(queryID[1], GL_TIMESTAMP);
// GL driver info and GLSL log
scr.text_init_pix(0.75);
glColor4f(1.0,1.0,1.0,0.9);
scr.text(glGetAnsiString(GL_VENDOR));
scr.text(glGetAnsiString(GL_RENDERER));
scr.text("OpenGL ver: "+glGetAnsiString(GL_VERSION));
if (_multi) scr.text("Multi pass");
else scr.text("Single pass");
if (shd.log.Length()!=41)
for (int i=1;i<=shd.log.Length();) scr.text(str_load_lin(shd.log,i,true));
scr.text_exit();
scr.exe();
scr.rfs();
// wait until the results are available
int e;
unsigned __int64 t0,t1;
for (e=0;!e;) glGetQueryObjectiv(queryID[0],GL_QUERY_RESULT_AVAILABLE,&e);
for (e=0;!e;) glGetQueryObjectiv(queryID[1],GL_QUERY_RESULT_AVAILABLE,&e);
glGetQueryObjectui64v(queryID[0], GL_QUERY_RESULT, &t0);
glGetQueryObjectui64v(queryID[1], GL_QUERY_RESULT, &t1);
Form1->Caption=AnsiString().sprintf("dt: %f ms p0:%.3fx%.3f zoom: %.1lf N:%i<<%i\n",(t1-t0)/1000000.0,mx,my,zoom,N,sh);
//---------------------------------------------------------------------------
__fastcall TForm1::TForm1(TComponent* Owner):TForm(Owner)
scr.init(this);
shd.set_source_file("","","","Mandelbrot_set.glsl_vert","Mandelbrot_set.glsl_frag");
glGenQueries(2, queryID);
// nice spirals
_multi=1;
zoom=300.0;
mx = 0.268;
my =-0.102;
//---------------------------------------------------------------------------
void __fastcall TForm1::FormDestroy(TObject *Sender)
scr.exit();
//---------------------------------------------------------------------------
void __fastcall TForm1::FormResize(TObject *Sender)
scr.resize();
xs=ClientWidth;
ys=ClientHeight;
txr.resize(xs,ys);
gl_draw();
//---------------------------------------------------------------------------
void __fastcall TForm1::FormPaint(TObject *Sender)
gl_draw();
//---------------------------------------------------------------------------
void __fastcall TForm1::FormMouseMove(TObject *Sender, TShiftState Shift, int X,int Y)
bool q0,q1;
mx1=1.0-divide(X+X,xs-1);
my1=divide(Y+Y,ys-1)-1.0;
sh1=Shift;
q0=sh0.Contains(ssLeft);
q1=sh1.Contains(ssLeft);
if (q1)
mx-=(mx1-mx0)/zoom;
my-=(my1-my0)/zoom;
mx0=mx1; my0=my1; sh0=sh1;
gl_draw();
//---------------------------------------------------------------------------
void __fastcall TForm1::FormMouseDown(TObject *Sender, TMouseButton Button,TShiftState Shift, int X, int Y)
FormMouseMove(Sender,Shift,X,Y);
//---------------------------------------------------------------------------
void __fastcall TForm1::FormMouseUp(TObject *Sender, TMouseButton Button,TShiftState Shift, int X, int Y)
FormMouseMove(Sender,Shift,X,Y);
//---------------------------------------------------------------------------
void __fastcall TForm1::FormMouseWheel(TObject *Sender, TShiftState Shift, int WheelDelta, TPoint &MousePos, bool &Handled)
if (WheelDelta>0) zoom*=1.2;
if (WheelDelta<0) zoom/=1.2;
Handled=true;
gl_draw();
//---------------------------------------------------------------------------
void __fastcall TForm1::FormKeyDown(TObject *Sender, WORD &Key, TShiftState Shift)
Caption=Key;
if (Key==32) _multi=!_multi; gl_draw(); // [Space]
if (Key==33) if (N<8192) N<<=1; gl_draw(); // [PgUp]
if (Key==34) if (N> 128) N>>=1; gl_draw(); // [PgDown]
//---------------------------------------------------------------------------
这是单程分数转义n=100*32：
这是单遍整数转义n=100：
正如您所见，对于相同的迭代次数，分数转义要好得多 (100)。
最后，只有 256 次迭代和约 300 倍缩放：
与单次通过：
关于修改的一些说明：
我将sh 小数位部分添加到计数器（定点）。所以现在最大计数是n<<sh 而不仅仅是n。我还添加了n0 常量，以降低转义小数部分的错误。该链接建议使用 2 次迭代，但我认为 1 次看起来更好（它还从对数方程中删除了 i+1 增量）。迭代循环没有改变，我只是在它之后添加相同的 n0 迭代，然后计算分数转义 mu 并将其转换为定点（因为我的着色器输出整数）。
仅在 CPU 端代码上更改了多通道。它只是重新索引使用过的索引，因此它们中没有孔，并使用可见光谱颜色重新着色。
这里演示：
Win32 C++/GL/GLSL Mandelbrot set demo 32 bit float (old) Win32 C++/GL/GLSL Mandelbrot set demo 64 bit float (new)【讨论】：
这个答案看起来很有希望，但我不太明白 Z(n) 是什么，请您详细介绍一下 Z(n) 是什么？ @MorganS42 正如我所写的，我还没有实现这个，所以我可能错了，但 Mandelbrot 集本身正在计算复域上的二次方程。如果您查看我的答案中的第一个链接，我在那里得到了我的 GLSL 实现。Z 是复杂域子结果，所以我假设在我的代码中它将是 Z=(x,y) 所以 |Z| = sqrt(x*x + y*y) 其中 x,y 值是在片段着色器中转义 for 循环后获取的......
是的，我确实看到它看起来非常好，但您认为您可以将公式应用于我提供的代码吗？
@MorganS42 您没有提供相关代码，所以不...更改将在mb(x, y, m) 函数中...您还调用mb(x, y, m) 两次？那将非常慢为什么不将结果存储到变量中呢？顺便说一句，通过搜索 fract_esc 宏，代码的更改很容易在我的代码中找到。我使用了 5 位小数的定点精度，因为我的着色器确实返回整数而不是浮点数 ...
@MorganS42 我添加了包含多通道代码的edit2（我设法让它与分数转义一起工作，它看起来很棒），我还发布了一个win32演示。新的被移植到double 并使用两个着色器来提高速度并添加了更多功能。你还是没有贴mb功能代码以上是关于找不到按照我的目标为 Mandelbrot 上色的方法的主要内容，如果未能解决你的问题，请参考以下文章