fixup bbh

tools/Status/ExecutableStatus/EvolveGhBinaryBlackHole.py

@@ -117,42 +117,42 @@ def render_dashboard(self, job: dict, input_file: dict):
         run_dir = Path(job["WorkDir"])
         reduction_files = list_reduction_files(job=job, input_file=input_file)
 
+        # Common horizon
+        with h5py.File(run_dir / reduction_files[-1], "r") as open_h5file:
+            ahc_subfile = open_h5file.get("ObservationAhC.dat")
+            if ahc_subfile is not None:
+                ahc_data = to_dataframe(ahc_subfile).iloc[0]
+                st.success(
+                    "Found a common horizon with mass"
+                    f" {ahc_data['ChristodoulouMass']:g} and spin"
+                    f" {ahc_data['DimensionlessSpinMagnitude']:g} at time"
+                    f" {ahc_data['Time']:g} M."
+                )
+
+        # Plot style
+        legend_layout = dict(
+            title=None,
+            orientation="h",
+            yanchor="bottom",
+            y=1,
+            xanchor="left",
+            x=0,
+        )
+
         # Constraints
         st.subheader("Constraints")
 
         def get_constraints_data(reductions_file):
             with h5py.File(reductions_file, "r") as open_h5file:
-                return pd.concat(
-                    [
-                        to_dataframe(
-                            open_h5file["ConstraintEnergy.dat"]
-                        ).set_index("Time")[["L2Norm(ConstraintEnergy)"]],
-                        to_dataframe(open_h5file["Norms.dat"]).set_index(
-                            "Time"
-                        )[
-                            [
-                                "L2Norm(PointwiseL2Norm(GaugeConstraint))",
-                                "L2Norm(PointwiseL2Norm(ThreeIndexConstraint))",
-                                "L2Norm(PointwiseL2Norm(FourIndexConstraint))",
-                            ]
-                        ],
-                    ],
-                    axis=1,
-                )
+                return to_dataframe(open_h5file["Norms.dat"]).set_index("Time")
 
         constraints = pd.concat(map(get_constraints_data, reduction_files))
         constraints.sort_index(inplace=True)
+        constraints = constraints[
+            [col for col in constraints.columns if "Constraint" in col]
+        ]
         fig = px.line(constraints.iloc[1:], log_y=True)
-        fig.update_layout(
-            legend=dict(
-                title=None,
-                orientation="h",
-                yanchor="bottom",
-                y=1,
-                xanchor="left",
-                x=0,
-            )
-        )
+        fig.update_layout(xaxis_title="Time [M]", legend=legend_layout)
         fig.update_yaxes(exponentformat="e", title=None)
         st.plotly_chart(fig)
 
@@ -161,15 +161,17 @@ def get_constraints_data(reductions_file):
 
         # Trajectories
         st.subheader("Trajectories")
-        st.pyplot(
-            plot_trajectory(
-                *import_A_and_B(
-                    reduction_files,
-                    "ApparentHorizons/ControlSystemAhA_Centers.dat",
-                    "ApparentHorizons/ControlSystemAhB_Centers.dat",
-                )
-            )
+        traj_A, traj_B = import_A_and_B(reduction_files)
+        st.pyplot(plot_trajectory(traj_A, traj_B))
+        coord_separation = np.linalg.norm(traj_A[:, 1:] - traj_B[:, 1:], axis=1)
+        fig = px.line(
+            x=traj_A[:, 0],
+            y=coord_separation,
+            labels={"y": "Coordinate separation [M]"},
         )
+        fig.update_layout(xaxis_title="Time [M]", legend=legend_layout)
+        fig.update_yaxes(title=None)
+        st.plotly_chart(fig)
 
         # Grid
         st.subheader("Grid")
@@ -212,26 +214,28 @@ def render_control_systems():
         render_control_systems()
 
         # Eccentricity
-        # (can be enabled once EccentricityControl is available)
-        # st.subheader("Eccentricity")
-
-        # @st.experimental_fragment
-        # def render_eccentricity():
-        #     if st.checkbox("Show eccentricity"):
-        #       ecc_control_result = coordinate_separation_eccentricity_control(
-        #             reduction_files[0],
-        #             "ApparentHorizons/ControlSystemAhA_Centers.dat",
-        #             "ApparentHorizons/ControlSystemAhB_Centers.dat",
-        #             tmin=60,
-        #             tmax=1000,
-        #             angular_velocity_from_xcts=None,
-        #             expansion_from_xcts=None,
-        #             output="temp.pdf",
-        #         )["H4"]["fit result"]
-        #         st.pyplot(plt.gcf())
-        #         st.metric(
-        #             "Eccentricity", f"{ecc_control_result['eccentricity']:e}"
-        #         )
-        #         st.write(ecc_control_result["xcts updates"])
-
-        # render_eccentricity()
+        st.subheader("Eccentricity")
+
+        @st.experimental_fragment
+        def render_eccentricity():
+            if st.checkbox("Show eccentricity"):
+                col_tmin, col_tmax = st.columns(2)
+                fig = plt.figure(figsize=(10, 8), layout="tight")
+                ecc_control_result = coordinate_separation_eccentricity_control(
+                    reduction_files[0],
+                    "ApparentHorizons/ControlSystemAhA_Centers.dat",
+                    "ApparentHorizons/ControlSystemAhB_Centers.dat",
+                    tmin=col_tmin.number_input("tmin", value=600, min_value=0),
+                    tmax=col_tmax.number_input("tmax", value=2000, min_value=0),
+                    angular_velocity_from_xcts=None,
+                    expansion_from_xcts=None,
+                    fig=fig,
+                )["H4"]["fit result"]
+                st.pyplot(fig)
+                st.metric(
+                    "Eccentricity",
+                    f"{ecc_control_result['eccentricity']:e}",
+                )
+                st.write(ecc_control_result)
+
+        render_eccentricity()