new Conformal Square projection (with transition and shader and Euclidean form)
This commit is contained in:
Zeno Rogue
parent
02f0b1e714
commit
1b7f4b869e
|
|
@ -1014,9 +1014,9 @@ enum eModel : int {
|
|||
mdWerner, mdAitoff, mdHammer, mdLoximuthal, mdMiller, mdGallStereographic, mdWinkelTripel,
|
||||
// 39..48
|
||||
mdPoorMan, mdPanini, mdRetroCraig, mdRetroLittrow, mdRetroHammer, mdThreePoint, mdLiePerspective, mdLieOrthogonal, mdRelPerspective, mdRelOrthogonal,
|
||||
// 49
|
||||
mdHorocyclicEqa,
|
||||
// 50..
|
||||
// 49..50
|
||||
mdHorocyclicEqa, mdConformalSquare,
|
||||
// 51..
|
||||
mdGUARD, mdPixel, mdHyperboloidFlat, mdPolynomial, mdManual
|
||||
};
|
||||
#endif
|
||||
|
|
@ -1078,6 +1078,7 @@ EX vector<modelinfo> mdinf = {
|
|||
{X3("relativistic perspective"), mf::euc_boring | mf::perspective},
|
||||
{X3("relativistic orthogonal"), mf::euc_boring},
|
||||
{X3("horocyclic equal-area"), mf::euc_boring | mf::equiarea | mf::orientation | mf::horocyclic},
|
||||
{X3("conformal square"), mf::orientation | mf::broken | mf::transition},
|
||||
{X3("guard"), mf::technical},
|
||||
{X3("pixel"), mf::technical},
|
||||
{X3("hypflat"), mf::technical},
|
||||
|
|
|
|||
56
drawing.cpp
56
drawing.cpp
|
|
@ -1622,10 +1622,14 @@ EX namespace ods {
|
|||
bool broken_projection(dqi_poly& p0) {
|
||||
int broken_coord = models::get_broken_coord(pmodel);
|
||||
static bool in_broken = false;
|
||||
bool both_broken = pmodel == mdConformalSquare;
|
||||
if(broken_coord && !in_broken) {
|
||||
|
||||
int zcoord = broken_coord;
|
||||
int ycoord = 3 - zcoord;
|
||||
int xcoord = 0;
|
||||
|
||||
zcoord = 2;
|
||||
|
||||
vector<hyperpoint> all;
|
||||
for(int i=0; i<p0.cnt; i++)
|
||||
|
|
@ -1634,9 +1638,17 @@ bool broken_projection(dqi_poly& p0) {
|
|||
int last_fail;
|
||||
|
||||
for(auto& h: all) models::apply_orientation(h[0], h[1]);
|
||||
|
||||
|
||||
auto break_in_xz = [&] (hyperpoint a, hyperpoint b, int xcoord, int zcoord) {
|
||||
return a[xcoord] * b[xcoord] <= 0 && (a[xcoord] * b[zcoord] - b[xcoord] * a[zcoord]) * (a[xcoord] - b[xcoord]) < 0;
|
||||
};
|
||||
|
||||
auto break_in = [&] (hyperpoint a, hyperpoint b) {
|
||||
return a[0] * b[0] <= 0 && (a[0] * b[zcoord] - b[0] * a[zcoord]) * (a[0] - b[0]) < 0;
|
||||
if(both_broken) {
|
||||
for(int xc=0; xc<2; xc++) {if(break_in_xz(a, b, xc, zcoord)) { xcoord = xc; ycoord = 1-xc; return true; } }
|
||||
return false;
|
||||
}
|
||||
return break_in_xz(a, b, xcoord, zcoord);
|
||||
};
|
||||
|
||||
for(int i=0; i<p0.cnt-1; i++)
|
||||
|
|
@ -1655,7 +1667,12 @@ bool broken_projection(dqi_poly& p0) {
|
|||
if(fail) {
|
||||
if(p0.tinf) return true;
|
||||
dynamicval<bool> ib(in_broken, true);
|
||||
ld part = ilerp(all[last_fail][0], all[last_fail+1][0], 0);
|
||||
|
||||
ld part = ilerp(all[last_fail][xcoord], all[last_fail+1][xcoord], 0);
|
||||
if(both_broken && all[last_fail][ycoord] * all[last_fail+1][ycoord] < 0) {
|
||||
ld part2 = ilerp(all[last_fail][ycoord], all[last_fail+1][ycoord], 0);
|
||||
if(part2 > part) part = part2, swap(xcoord, ycoord);
|
||||
}
|
||||
hyperpoint initial = normalize(lerp(all[last_fail], all[last_fail+1], 1 - (1-part) * .99));
|
||||
bool have_initial = true;
|
||||
v.push_back(glhr::pointtogl(initial));
|
||||
|
|
@ -1663,19 +1680,46 @@ bool broken_projection(dqi_poly& p0) {
|
|||
int at = last_fail;
|
||||
do {
|
||||
v.push_back(glhr::pointtogl(all[at]));
|
||||
if(at == p0.cnt-1 && all[at] != all[0]) {
|
||||
if(at == p0.cnt-1 && sqhypot_d(2, all[at] - all[0]) > 1e-6) {
|
||||
p.cnt = isize(v); p.draw(); v.clear(); at = 0;
|
||||
have_initial = false;
|
||||
}
|
||||
int next = at+1;
|
||||
if(next == p0.cnt) next = 0;
|
||||
if(break_in(all[at], all[next])) {
|
||||
ld part = ilerp(all[at][0], all[next][0], 0);
|
||||
ld part = ilerp(all[at][xcoord], all[next][xcoord], 0);
|
||||
if(both_broken && all[at][ycoord] * all[next][ycoord] < 0) {
|
||||
ld part2 = ilerp(all[at][ycoord], all[next][ycoord], 0);
|
||||
if(part2 < part) part = part2, swap(xcoord, ycoord);
|
||||
}
|
||||
|
||||
hyperpoint final = normalize(lerp(all[at], all[next], part * .99));
|
||||
v.push_back(glhr::pointtogl(final));
|
||||
if(have_initial) {
|
||||
int max = 4 << vid.linequality;
|
||||
if(final[0] * initial[0] > 0) {
|
||||
if(both_broken) {
|
||||
auto square_close_corner = [] (hyperpoint h) {
|
||||
hyperpoint end = -C0;
|
||||
if(abs(h[0]) > abs(h[1]))
|
||||
end[0] = 0.01 * signum(h[0]), end[1] = 0.001 * signum(h[1]);
|
||||
else
|
||||
end[1] = 0.01 * signum(h[1]), end[0] = 0.001 * signum(h[0]);
|
||||
return normalize(end);
|
||||
};
|
||||
hyperpoint endf = square_close_corner(final);
|
||||
hyperpoint endi = square_close_corner(initial);
|
||||
if(endf != endi) {
|
||||
for(int i=1; i<=max; i++)
|
||||
v.push_back(glhr::pointtogl(lerp(final, endf, i * 1. / max)));
|
||||
for(int i=0; i<=max; i++)
|
||||
v.push_back(glhr::pointtogl(lerp(endi, initial, i * 1. / max)));
|
||||
}
|
||||
else {
|
||||
for(int i=1; i<=max; i++)
|
||||
v.push_back(glhr::pointtogl(lerp(final, initial, i * 1. / max)));
|
||||
}
|
||||
}
|
||||
else if(final[xcoord] * initial[xcoord] > 0) {
|
||||
for(int i=1; i<=max; i++)
|
||||
v.push_back(glhr::pointtogl(lerp(final, initial, i * 1. / max)));
|
||||
}
|
||||
|
|
|
|||
85
hypgraph.cpp
85
hypgraph.cpp
|
|
@ -528,6 +528,81 @@ EX void make_axial(hyperpoint H, hyperpoint& ret, const hr::function<ld(hyperpoi
|
|||
models::apply_orientation_yz(ret[2], ret[1]);
|
||||
}
|
||||
|
||||
// according to https://github.com/cspersonal/peirce-quincuncial-projection/blob/master/peirceQuincuncialProjection.R
|
||||
|
||||
ld ellRF(ld x, ld y, ld z) {
|
||||
ld delx = 1, dely = 1, delz = 1;
|
||||
const ld eps = 0.0025;
|
||||
ld mean;
|
||||
while(abs(delx) > eps || abs(dely) > eps || abs(delz) > eps) {
|
||||
ld sx = sqrt(x);
|
||||
ld sy = sqrt(y);
|
||||
ld sz = sqrt(z);
|
||||
ld len = sx * (sy+sz) + sy * sz;
|
||||
x = .25 * (x+len);
|
||||
y = .25 * (y+len);
|
||||
z = .25 * (z+len);
|
||||
mean = (x+y+z)/3;
|
||||
delx = (mean-x) / mean;
|
||||
dely = (mean-y) / mean;
|
||||
delz = (mean-z) / mean;
|
||||
}
|
||||
ld e2 = delx * dely - delz * delz;
|
||||
ld e3 = delx * dely * delz;
|
||||
return ((1.0 + (e2 / 24.0 - 0.1 - 3.0 * e3 / 44.0) * e2+ e3 / 14) / sqrt(mean));
|
||||
}
|
||||
|
||||
ld ellFaux(ld cos_phi, ld sin_phi, ld k) {
|
||||
ld x = cos_phi * cos_phi;
|
||||
ld y = 1 - k * k * sin_phi * sin_phi;
|
||||
return sin_phi * ellRF(x, y, 1);
|
||||
}
|
||||
|
||||
ld sqrt_clamp(ld x) { if(x<0) return 0; return sqrt(x); }
|
||||
|
||||
hyperpoint to_square(hyperpoint H) {
|
||||
|
||||
ld d = hypot_d(2, H);
|
||||
ld x = d / (H[2] + 1);
|
||||
x *= pconf.model_transition;
|
||||
|
||||
ld cos_phiosqrt2 = sqrt(2) / (x + 1/x);
|
||||
ld cos_lambda = -H[1] / d;
|
||||
ld sin_lambda = H[0] / d;
|
||||
ld cos_a = cos_phiosqrt2 * (sin_lambda + cos_lambda);
|
||||
ld cos_b = cos_phiosqrt2 * (sin_lambda - cos_lambda);
|
||||
ld sin_a = sqrt(1 - cos_a * cos_a);
|
||||
ld sin_b = sqrt(1 - cos_b * cos_b);
|
||||
ld cos_a_cos_b = cos_a * cos_b;
|
||||
ld sin_a_sin_b = sin_a * sin_b;
|
||||
ld sin2_m = 1.0 + cos_a_cos_b - sin_a_sin_b;
|
||||
ld sin2_n = 1.0 - cos_a_cos_b - sin_a_sin_b;
|
||||
ld sin_m = sqrt_clamp(sin2_m);
|
||||
ld cos_m = sqrt_clamp(1 - sin2_m);
|
||||
if(sin_lambda < 0) sin_m = -sin_m;
|
||||
ld sin_n = sqrt_clamp(sin2_n);
|
||||
ld cos_n = sqrt_clamp(1.0 - sin2_n);
|
||||
if(cos_lambda > 0.0) sin_n = -sin_n;
|
||||
|
||||
hyperpoint res;
|
||||
ld divby = 0.53935260118837935472;
|
||||
res[0] = ellFaux(cos_m,sin_m,sqrt(2)/2.) * divby;
|
||||
res[1] = ellFaux(cos_n,sin_n,sqrt(2)/2.) * divby;
|
||||
res[2] = 0; res[3] = 1;
|
||||
|
||||
if(x > 1) {
|
||||
if(abs(res[0]) > abs(res[1])) {
|
||||
if(res[0] > 0) res[0] = 2 - res[0]; else res[0] = -2 - res[0];
|
||||
}
|
||||
else {
|
||||
if(res[1] > 0) res[1] = 2 - res[1]; else res[1] = -2 - res[1];
|
||||
}
|
||||
}
|
||||
|
||||
res /= pconf.model_transition;
|
||||
return res;
|
||||
}
|
||||
|
||||
EX void apply_other_model(shiftpoint H_orig, hyperpoint& ret, eModel md) {
|
||||
|
||||
hyperpoint H = H_orig.h;
|
||||
|
|
@ -810,6 +885,16 @@ EX void apply_other_model(shiftpoint H_orig, hyperpoint& ret, eModel md) {
|
|||
break;
|
||||
}
|
||||
|
||||
case mdConformalSquare: {
|
||||
find_zlev(H);
|
||||
models::apply_orientation_yz(H[1], H[2]);
|
||||
models::apply_orientation(H[0], H[1]);
|
||||
ret = to_square(H);
|
||||
models::apply_orientation(ret[1], ret[0]);
|
||||
models::apply_orientation_yz(ret[2], ret[1]);
|
||||
break;
|
||||
}
|
||||
|
||||
case mdLieOrthogonal: {
|
||||
ret = lie_log_correct(H_orig, H);
|
||||
|
||||
|
|
|
|||
68
shaders.cpp
68
shaders.cpp
|
|
@ -276,10 +276,15 @@ shared_ptr<glhr::GLprogram> write_shader(flagtype shader_flags) {
|
|||
}
|
||||
else if(pmodel == mdDisk && GDIM == 3 && !spherespecial && !nonisotropic && !gproduct) {
|
||||
coordinator += "t /= (t[3] + uAlpha);\n";
|
||||
vsh += "uniform mediump float uAlpha;";
|
||||
vsh += "uniform mediump float uAlpha;\n";
|
||||
shader_flags |= SF_DIRECT | SF_BOX | SF_ZFOG;
|
||||
treset = true;
|
||||
}
|
||||
else if(pmodel == mdConformalSquare && pconf.model_transition == 1) {
|
||||
shader_flags |= SF_ORIENT | SF_DIRECT;
|
||||
coordinator += "t = uPP * t;", vsh += "uniform mediump mat4 uPP;";
|
||||
coordinator += "t = to_square(t);";
|
||||
}
|
||||
else if(pmodel == mdBand && hyperbolic) {
|
||||
shader_flags |= SF_BAND | SF_ORIENT | SF_BOX | SF_DIRECT;
|
||||
coordinator += "t = uPP * t;", vsh += "uniform mediump mat4 uPP;";
|
||||
|
|
@ -848,6 +853,67 @@ EX void add_if(string& shader, const string& seek, const string& function) {
|
|||
EX void add_fixed_functions(string& shader) {
|
||||
/* from the most complex to the simplest */
|
||||
|
||||
add_if(shader, "to_square",
|
||||
"mediump vec4 to_square(mediump vec4 h) {\n"
|
||||
"float d = length(h.xy);\n"
|
||||
"float x = d / (h.z + 1.);\n"
|
||||
|
||||
"float cos_phiosqrt2 = sqrt(2.) / (x + 1./x);\n"
|
||||
"float cos_lambda = -h.y / d;\n"
|
||||
"float sin_lambda = h.x / d;\n"
|
||||
"float cos_a = cos_phiosqrt2 * (sin_lambda + cos_lambda);\n"
|
||||
"float cos_b = cos_phiosqrt2 * (sin_lambda - cos_lambda);\n"
|
||||
"float sin_a = sqrt(1. - cos_a * cos_a);\n"
|
||||
"float sin_b = sqrt(1. - cos_b * cos_b);\n"
|
||||
"float cos_a_cos_b = cos_a * cos_b;\n"
|
||||
"float sin_a_sin_b = sin_a * sin_b;\n"
|
||||
"float sin2_m = 1.0 + cos_a_cos_b - sin_a_sin_b;\n"
|
||||
"float sin2_n = 1.0 - cos_a_cos_b - sin_a_sin_b;\n"
|
||||
"float sin_m = sqrt_clamp(sin2_m);\n"
|
||||
"float cos_m = sqrt_clamp(1. - sin2_m);\n"
|
||||
"if(sin_lambda < 0.) sin_m = -sin_m;\n"
|
||||
"float sin_n = sqrt_clamp(sin2_n);\n"
|
||||
"float cos_n = sqrt_clamp(1.0 - sin2_n);\n"
|
||||
"if(cos_lambda > 0.0) sin_n = -sin_n;\n"
|
||||
"#define divby 0.53935260118837935472\n"
|
||||
"vec4 res = vec4(ellFaux(cos_m,sin_m,sqrt(2.)/2.) * divby, ellFaux(cos_n,sin_n,sqrt(2.)/2.) * divby, 0, 1);\n"
|
||||
"if(x > 1.) {\n"
|
||||
" if(abs(res[0]) > abs(res[1])) {\n"
|
||||
" if(res[0] > 0.) res[0] = 2. - res[0]; else res[0] = -2. - res[0];\n"
|
||||
" }\n"
|
||||
" else {\n"
|
||||
" if(res[1] > 0.) res[1] = 2. - res[1]; else res[1] = -2. - res[1];\n"
|
||||
" }\n"
|
||||
" }\n"
|
||||
"return res;\n"
|
||||
"}\n");
|
||||
|
||||
add_if(shader, "sqrt_clamp", "mediump float sqrt_clamp(mediump float x) { return x >= 0. ? sqrt(x) : 0.; }\n");
|
||||
add_if(shader, "ellFaux", "mediump float ellFaux(mediump float cos_phi, mediump float sin_phi, mediump float k) {\n"
|
||||
"return sin_phi * ellRF(cos_phi * cos_phi, 1. - k * k * sin_phi * sin_phi, 1.);\n"
|
||||
"}\n");
|
||||
add_if(shader, "ellRF", "mediump float ellRF(mediump float x, mediump float y, mediump float z) {\n"
|
||||
"float delx = 1., dely = 1., delz = 1.;\n"
|
||||
"const float eps = 0.0025;\n"
|
||||
"float mean;\n"
|
||||
"while(abs(delx) > eps || abs(dely) > eps || abs(delz) > eps) {\n"
|
||||
" float sx = sqrt(x);\n"
|
||||
" float sy = sqrt(y);\n"
|
||||
" float sz = sqrt(z);\n"
|
||||
" float len = sx * (sy+sz) + sy * sz;\n"
|
||||
" float x = .25 * (x+len);\n"
|
||||
" float y = .25 * (y+len);\n"
|
||||
" float z = .25 * (z+len);\n"
|
||||
" mean = (x+y+z)/3.;\n"
|
||||
" delx = (mean-x) / mean;\n"
|
||||
" dely = (mean-y) / mean;\n"
|
||||
" delz = (mean-z) / mean;\n"
|
||||
" }\n"
|
||||
"float e2 = delx * dely - delz * delz;\n"
|
||||
"float e3 = delx * dely * delz;\n"
|
||||
"return ((1.0 + (e2 / 24.0 - 0.1 - 3.0 * e3 / 44.0) * e2+ e3 / 14.) / sqrt(mean));\n"
|
||||
"}\n");
|
||||
|
||||
add_if(shader, "tanh", "mediump float tanh(mediump float x) { return sinh(x) / cosh(x); }\n");
|
||||
add_if(shader, "sinh", "mediump float sinh(mediump float x) { return (exp(x) - exp(-x)) / 2.0; }\n");
|
||||
add_if(shader, "asin_clamp", "mediump float asin_clamp(mediump float x) { return x > 1. ? PI/2. : x < -1. ? -PI/2. : asin(x); }\n");
|
||||
|
|
|
|||
Loading…
Reference in New Issue