new embedding: Euclidean cylinder

3d-models.cpp

@@ -35,6 +35,7 @@ hyperpoint get_center(const vector<hyperpoint>& vh) {
   hyperpoint h = Hypc;
   for(auto h1: vh) h = h + h1;
   if(geom3::euc_in_product()) return h / isize(vh);
+  if(geom3::euc_cylinder()) h /= isize(vh);
   return normalize_flat(h);
   }
 

geometry.cpp

@@ -599,6 +599,7 @@ void geometry_information::prepare_lta() {
   if(embedded_plane) {
     if(geom3::euc_vertical()) lta = cspin90(2, 1) * lta;
     if(geom3::hyp_in_solnih()) lta = cspin90(0, 1) * cspin90(1, 2) * cspin90(0, 1) * lta;
+    if(geom3::euc_cylinder()) lta = cspin90(0, 1) * lta;
     }
   logical_scaled_to_intemediate = lta;
   if(geom3::euc_in_noniso()) {
@@ -1087,7 +1088,11 @@ EX namespace geom3 {
       HIGH2 = lev_to_factor(3 * wh);
       SKY = LOWSKY - sgn * 5;
 
-      if(geom3::mgclass() == gcSphere && geom3::ggclass() != gcSphere) {
+      /* in spherical/cylindrical case, make sure that the high stuff does not go through the center */
+      bool depth_limit = geom3::mgclass() == gcSphere && geom3::ggclass() != gcSphere;
+      depth_limit |= euc_cylinder();
+
+      if(depth_limit) {
         ld max_high = lerp(-FLOOR, -1, 0.8);
         ld max_high2 = lerp(-FLOOR, -1, 0.9);
         if(HIGH < max_high) HIGH = max_high;
@@ -1123,7 +1128,8 @@ EX namespace geom3 {
     seCliffordTorus,
     seProductH,
     seProductS,
-    seSL2
+    seSL2,
+    seCylinder
     };
   #endif
 
@@ -1140,7 +1146,8 @@ EX namespace geom3 {
     {"Clifford Torus",              "Embed Euclidean rectangular torus into S3."},
     {"hyperbolic product", "Embed Euclidean or hyperbolic plane in the H2xR product space."},
     {"spherical product", "Embed Euclidean cylinder or spherical plane in the H2xR product space."},
-    {"SL(2,R)",           "Embed Euclidean plane in twisted product geometry."}
+    {"SL(2,R)",           "Embed Euclidean plane in twisted product geometry."},
+    {"cylinder",          "Embed Euclidean cylinder in Euclidean space."},
     };
 
   EX eSpatialEmbedding spatial_embedding = seDefault;
@@ -1196,9 +1203,14 @@ EX namespace geom3 {
     }
 
   EX bool euc_in_noniso() {
+    if(spatial_embedding == seCylinder) return mgclass() == gcEuclid;
     return among(ggclass(), gcNil, gcSol, gcNIH, gcSolN, gcSphere, gcProduct, gcSL2) && mgclass() == gcEuclid;
     }
 
+  EX bool euc_cylinder() {
+    return spatial_embedding == seCylinder && mgclass() == gcEuclid;
+    }
+
   EX bool sph_in_euc() {
     return ggclass() == gcEuclid && mgclass() == gcSphere;
     }
@@ -1220,7 +1232,7 @@ EX namespace geom3 {
     }
 
   EX bool same_in_same() {
-    return mgclass() == ggclass();
+    return mgclass() == ggclass() && !among(spatial_embedding, seCylinder);
     }
 
   EX bool flipped;
@@ -1450,7 +1462,8 @@ EX void switch_always3() {
           }
         }
       if(spatial_embedding == seCliffordTorus) configure_clifford_torus();
-      if(spatial_embedding == seProductS) configure_product_cylinder();
+      if(spatial_embedding == seProductS) configure_cylinder();
+      if(spatial_embedding == seCylinder) configure_cylinder();
       }
     else {
       vid.always3 = false;
@@ -1493,7 +1506,7 @@ EX void switch_always3() {
     vid.eye = vid.wall_height / 2 - vid.depth;
     }
 
-  EX void configure_product_cylinder() {
+  EX void configure_cylinder() {
     rug::clifford_torus ct;
     hyperpoint vec;
     if(sqhypot_d(2, ct.yh) > 1e-6) vec = ct.yh;

graph.cpp

@@ -785,6 +785,7 @@ EX shiftmatrix face_the_player(const shiftmatrix V) {
     }
   #endif
   if(embedded_plane && geom3::sph_in_low()) return shiftless(map_relative_push(unshift(V * zpush0(1))) * zpush(-1));
+  if(embedded_plane && geom3::euc_cylinder()) return shiftless(map_relative_push(unshift(V * zpush0(1))) * zpush(-1));
   return rgpushxto0(tC0(V));
   }
 
@@ -5104,6 +5105,9 @@ EX void make_actual_view() {
     else if(geom3::euc_in_sph()) {
       current_display->radar_transform = inverse(View);
       }
+    else if(geom3::euc_cylinder()) {
+      current_display->radar_transform = inverse(View);
+      }
     else {
       transmatrix T = actual_view_transform * View;
       transmatrix U = view_inverse(T);

hyperpoint.cpp

@@ -629,6 +629,12 @@ EX hyperpoint normalize_flat(hyperpoint h) {
     h1[2] = 0;
     return parabolic13(h1);
     }
+  if(geom3::euc_cylinder()) {
+    h /= h[3];
+    ld z = h[1] * h[1] + h[2] * h[2];
+    if(z > 0) h[1] /= z, h[2] /= z;
+    return h;
+    }
   if(geom3::sph_in_euc()) {
     ld z = hypot_d(3, h);
     if(z > 0) h[0] /= z, h[1] /= z, h[2] /= z;
@@ -883,6 +889,14 @@ EX hyperpoint orthogonal_move(const hyperpoint& h, ld z) {
     hf[3] = 1;
     return hf;
     }
+  if(geom3::euc_cylinder()) {
+    auto hf = h / h[3];
+    ld z0 = hypot(h[1], h[2]);
+    if(!z0) return hf;
+    ld f = ((z0 + z) / z0);
+    hf[1] *= f; hf[2] *= f;
+    return hf;
+    }
   if(geom3::hyp_in_solnih()) {
     return nisot::translate(h) * cpush0(0, z);
     }
@@ -946,6 +960,9 @@ EX ld get_logical_z(hyperpoint h) {
     auto h1 = h / exp(bz);
     return asin_auto(h1[1]);
     }
+  if(geom3::euc_cylinder()) {
+    return hypot(h[1], h[2]) - 1;
+    }
   if(gproduct)
     return log(h[2]);
   return asin_auto(h[2]) - (moved_center() ? 1 : 0);
@@ -999,6 +1016,13 @@ EX void swapmatrix(transmatrix& T) {
       for(int i=0; i<4; i++) T[i][3] = T[3][i] = i == 3;
       }
     }
+  else if(geom3::euc_cylinder()) {
+    if(!geom3::flipped) {
+      hyperpoint h1 = cgi.logical_to_intermediate * get_column(T, 2);
+      T = xpush(h1[0]) * cspin(1, 2, h1[1]);
+      return;
+      }
+    }
   else if(geom3::euc_in_nil()) {
     if(!geom3::flipped) {
       hyperpoint h1 = cgi.logical_to_intermediate * get_column(T, 2);
@@ -1054,6 +1078,14 @@ EX void swapmatrix(hyperpoint& h) {
   if(geom3::in_product()) return;
   if(geom3::sph_in_euc()) { h[3] = 1; return; }
   if(geom3::sph_in_hyp()) { h[0] *= sinh(1); h[1] *= sinh(1); h[2] *= sinh(1); h[3] = cosh(1); return; }
+  if(geom3::euc_cylinder()) {
+    hyperpoint h1 = cgi.logical_to_intermediate * h;
+    h[0] = h1[0];
+    h[1] = sin(h1[1]);
+    h[2] = cos(h1[1]);
+    h[3] = 1;
+    return;
+    }
   if(geom3::euc_in_nil()) { h = cgi.logical_to_intermediate * h; h[3] = 1; h[1] = 0; return; }
   if(geom3::euc_in_sl2()) {
     hyperpoint h1 = cgi.logical_to_intermediate * h; h1[1] = 0;
@@ -1660,6 +1692,7 @@ EX bool moved_center() {
   if(geom3::sph_in_euc()) return true;
   if(geom3::sph_in_hyp()) return true;
   if(geom3::euc_in_sph()) return true;
+  if(geom3::euc_cylinder()) return true;
   return false;
   }
 
@@ -1668,6 +1701,7 @@ EX hyperpoint tile_center() {
   if(geom3::sph_in_euc()) return C02 + C03;
   if(geom3::euc_in_sph()) return zpush0(1);
   if(geom3::sph_in_hyp()) return zpush0(1);
+  if(geom3::euc_cylinder()) return zpush0(1);
   return C0;
   }
 

hypgraph.cpp

@@ -3398,6 +3398,12 @@ EX transmatrix map_relative_push(hyperpoint h) {
     auto h1 = esl2_ati(h);
     return esl2_zpush(h1[2]) * xpush(h1[0]) * ypush(h1[1]);
     }
+  if(geom3::euc_cylinder()) {
+    ld z0 = hypot(h[1], h[2]);
+    if(!z0) return Id;
+    transmatrix T = xpush(h[0]) * cspin(1, 2, atan2(h[1], h[2])) * zpush(z0);
+    return T;
+    }
   if(geom3::euc_in_sph()) {
     ld tx = hypot(h[0], h[2]);
     ld ty = hypot(h[1], h[3]);

radar.cpp

@@ -63,6 +63,15 @@ pair<bool, hyperpoint> makeradar(shiftpoint h) {
       if(d > vid.radarrange) return {false, h1};
       if(d) h1 = h1 / (vid.radarrange + cgi.scalefactor/4);
       }
+    else if(geom3::euc_cylinder()) {
+      h1[0] = h.h[0];
+      h1[1] = atan2(h.h[1], h.h[2]);
+      h1[2] = 0;
+      h1 = cgi.intermediate_to_logical * h1;
+      d = hypot_d(2, h1);
+      if(d > vid.radarrange) return {false, h1};
+      if(d) h1 = h1 / (vid.radarrange + cgi.scalefactor/4);
+      }
     else if(geom3::euc_in_sl2()) {
       h1 = cgi.intermediate_to_logical * esl2_ati(unshift(h)); h1[1] = -h1[1];
       d = hypot_d(2, h1);

sky.cpp

@@ -43,7 +43,7 @@ EX struct dqi_sky *sky;
 
 EX void prepare_sky() {
   sky = NULL;
-  if(euclid && !geom3::sph_in_euc()) {
+  if(euclid && !geom3::sph_in_euc() && !geom3::euc_cylinder()) {
     if(WDIM == 3 || GDIM == 2) return;
     if(no_wall_rendering) return;
     if(!draw_sky) return;
@@ -76,6 +76,7 @@ void compute_skyvertices(const vector<sky_item>& sky) {
   if(geom3::hyp_in_solnih()) return;
   if(geom3::euc_in_product()) return;
   if(geom3::euc_in_sl2()) return;
+  if(geom3::euc_cylinder()) return;
 
   int sk = get_skybrightness();