Add direct cutting-plane control to GLVis scripts

lib/script_controller.cpp

@@ -20,6 +20,8 @@
 #include <array>
 #include <algorithm>
 #include <thread>
+#include <cmath>
+#include <cctype>
 
 using namespace std;
 using namespace mfem;
@@ -63,6 +65,7 @@ enum class Command
    PlotCaption,
    PointLine,
    Headless,
+   CuttingPlane,
    //----------
    Max
 };
@@ -128,6 +131,7 @@ ScriptCommands::ScriptCommands()
    (*this)[Command::PlotCaption]          = {"plot_caption", "'<caption>'", "Set the plot caption."};
    (*this)[Command::PointLine]            = {"pointline", "<num_points> <x y z>...", "Set point line overlay coordinates."};
    (*this)[Command::Headless]             = {"headless", "", "Change the session to headless."};
+   (*this)[Command::CuttingPlane]         = {"cutting_plane", "<phi> <theta> <translation> <kind:optional, default 1> <alg:optional>", "Set the cutting plane orientation (degrees), translation, kind, and algorithm."};
 }
 
 int ScriptController::ScriptReadSolution(istream &scr, DataState &state)
@@ -859,6 +863,49 @@ bool ScriptController::ExecuteScriptCommand()
          case Command::Headless:
             cout << "The session cannot become headless after initialization" << endl;
             break;
+         case Command::CuttingPlane:
+         {
+            double phi_deg, theta_deg, translation;
+            scr >> phi_deg >> theta_deg >> translation;
+
+            int kind = 1, algo = -1;
+            scr >> ws;
+            int ch = scr.peek();
+            if (isdigit(ch) || ch == '-')
+            {
+               scr >> kind;
+               scr >> ws;
+               ch = scr.peek();
+               if (isdigit(ch) || ch == '-')
+               {
+                  scr >> algo;
+               }
+            }
+
+            if (kind < -1 || kind > 2)
+            {
+               cerr << "Script: cutting_plane: invalid kind " << kind
+                    << " (expected -1..2)" << endl;
+               kind = 1;
+            }
+            if (algo != -1 && (algo < 0 || algo > 1))
+            {
+               cerr << "Script: cutting_plane: invalid alg " << algo
+                    << " (expected 0 or 1)" << endl;
+               algo = -1;
+            }
+
+            cout << "Script: cutting_plane: " << phi_deg << ' ' << theta_deg
+                 << ' ' << translation << ' ' << kind;
+            if (algo != -1) { cout << ' ' << algo; }
+            cout << endl;
+
+            win.vs->SetCuttingPlane(phi_deg * M_PI / 180.0,
+                                    theta_deg * M_PI / 180.0,
+                                    translation, kind, algo);
+            MyExpose();
+         }
+         break;
          case Command::Max: //dummy
             break;
       }

lib/vsdata.cpp

@@ -1986,6 +1986,7 @@ Plane::Plane(const double (&eqn_)[4], const VisualizationScene::Box &bb)
    x0 = (x[0]+x[1])/2.0;
    y0 = (y[0]+y[1])/2.0;
    z0 = (z[0]+z[1])/2.0;
+   cx = x0; cy = y0; cz = z0;
 
    phi_step = M_PI / 36;
    theta_step = M_PI / 36;
@@ -2050,3 +2051,19 @@ void Plane::DecreaseDistance()
    eqn[3] -= rho_step;
    CartesianToSpherical();
 }
+
+void Plane::SetPlane(double phi_, double theta_, double translation)
+{
+   phi = phi_;
+   theta = theta_;
+   rho = 1.0;
+
+   double nx = cos(phi) * cos(theta);
+   double ny = cos(phi) * sin(theta);
+   double nz = sin(phi);
+   x0 = cx + translation * nx;
+   y0 = cy + translation * ny;
+   z0 = cz + translation * nz;
+
+   SphericalToCartesian();
+}

lib/vsdata.hpp

@@ -23,6 +23,7 @@ class Plane
    double eqn[4];
    double phi, theta, rho;
    double x0,y0,z0;
+   double cx, cy, cz; // fixed mesh bounding-box center
    void CartesianToSpherical();
    void SphericalToCartesian();
 
@@ -44,6 +45,7 @@ class Plane
    void DecreaseTheta();
    void IncreaseDistance();
    void DecreaseDistance();
+   void SetPlane(double phi_, double theta_, double translation);
 };
 
 
@@ -263,6 +265,12 @@ class VisualizationSceneScalarData : public VisualizationScene
    virtual void AutoRefine() = 0;
    virtual void ToggleAttributes(mfem::Array<int> &attr_list) = 0;
 
+   // Set the cutting-plane orientation (radians), translation, kind, and
+   // algorithm in one shot. No-op by default; only meaningful for 3D scenes.
+   // kind == -1 / algo == -1 mean "leave unchanged".
+   virtual void SetCuttingPlane(double phi, double theta, double translation,
+                                int kind, int algo) { }
+
    virtual void PrintState();
 
    mfem::Mesh *GetMesh() { return mesh; }

lib/vssolution3d.cpp

@@ -1192,6 +1192,42 @@ void VisualizationSceneSolution3d::ToggleCPAlgorithm()
    }
 }
 
+void VisualizationSceneSolution3d::SetCuttingPlane(double phi, double theta,
+                                                   double translation,
+                                                   int kind, int algo)
+{
+   const int old_cplane = cplane;
+
+   CuttingPlane->SetPlane(phi, theta, translation);
+   FindNodePos();
+
+   if (kind != -1 && kind != cplane)
+   {
+      if (cplane == 2 && cp_drawmesh == 3)
+      {
+         cp_drawmesh = 2;
+      }
+      cplane = kind;
+   }
+
+   if (algo != -1)
+   {
+      cp_algo = algo;
+   }
+
+   CPPrepare();
+   if (old_cplane == 2 || cplane == 2)
+   {
+      Prepare();
+      PrepareLines();
+
+      if (cplane != old_cplane)
+      {
+         PrepareOrderingCurve();
+      }
+   }
+}
+
 void VisualizationSceneSolution3d::MoveLevelSurf(int move)
 {
    drawlsurf += move;

lib/vssolution3d.hpp

@@ -218,6 +218,8 @@ class VisualizationSceneSolution3d : public VisualizationSceneScalarData
    void ToggleCPDrawElems();
    void ToggleCPDrawMesh();
    void ToggleCPAlgorithm();
+   void SetCuttingPlane(double phi, double theta, double translation,
+                        int kind, int algo) override;
    void MoveLevelSurf(int);
    void NumberOfLevelSurf(int);
    void EventUpdateColors() override;