Update a script to check 1d senstivity kernel and 1d model

.gitignore

@@ -54,6 +54,7 @@ __pycache__/
 
 # Output files produced by the program
 OUTPUT_FILES/
+figures/
 
 tests/*/*.png
 tests/*/*.h5

examples/00_checkerboard_iso/run_this_example.sh

@@ -7,7 +7,7 @@ cp src_rec_file_ph.csv OUTPUT_FILES/src_rec_file_forward_RL_PH.csv
 # cp src_rec_file_ph.csv OUTPUT_FILES/src_rec_file_forward_LV_PH.csv
 
 # create 2x3x2 checkers and forward simulate surface traveltimes
-mpirun -np $NPROC ../../bin/surfatt_cb_fwd -i input_params.yml -n 2/3/2 -m 0.2 -p 0.08
+mpirun -np $NPROC ../../bin/SURFATT_cb_fwd -i input_params.yml -n 2/3/2 -m 0.2 -p 0.08
 
 # inversion 
-mpirun -np $NPROC ../../bin/surfatt_tomo -i input_params.yml
+mpirun -np $NPROC ../../bin/SURFATT_tomo -i input_params.yml

include/h5io.h

@@ -94,6 +94,9 @@ class H5IO {
     std::vector<T> read_vector(const std::string &name) const {
         H5::DataSet   ds = file_.openDataSet(name);
         H5::DataSpace sp = ds.getSpace();
+        if (sp.getSimpleExtentNdims() != 1)
+            throw std::runtime_error(
+                "H5IO::read_vector: dataset '" + name + "' is not 1-D");
         hsize_t n = 0;
         sp.getSimpleExtentDims(&n, nullptr);
         std::vector<T> v(n);

include/inversion1d.h

@@ -7,6 +7,8 @@
 #include <Eigen/Core>
 #include <Eigen/Dense>
 
+namespace surfker { struct DepthKernel1D; }
+
 class Inversion1D {
 public:
     explicit Inversion1D(WaveType wavetype);
@@ -26,4 +28,23 @@ class Inversion1D {
 
     int niter;
     std::vector<real_t> misfits;
+
+    // Validate depth grid, initial Vs and the active src_rec tables before
+    // the first forward call; logs every problem found and aborts if any.
+    void check_inputs(const Eigen::VectorX<real_t>& zarr) const;
+
+    // Abort with a diagnostic if disper() failed (returned 0 / non-finite).
+    void check_pred_vel(const Eigen::VectorX<real_t>& pred_vel,
+                        const Eigen::VectorX<real_t>& periods,
+                        const Eigen::VectorX<real_t>& zarr,
+                        int iter, WaveType wt, SurfType tp) const;
+
+    // Abort with a diagnostic if the depth kernels contain inf/NaN.
+    void check_kernels(const surfker::DepthKernel1D& kernels,
+                       const Eigen::VectorX<real_t>& periods,
+                       const Eigen::VectorX<real_t>& zarr,
+                       int iter, WaveType wt, SurfType tp) const;
+
+    // Log the current Vs profile and its first non-physical node, if any.
+    void log_vs_profile(const Eigen::VectorX<real_t>& zarr, int n_update) const;
 };

include/src_rec.h

@@ -72,7 +72,13 @@ class SrcRec {
     // Rank 0 reads the file; each field gets its own per-node MPI shared-memory
     // window via Parallel::alloc_shared (1-D for numeric, 2-D for strings).
     // Must call release_shm() before MPI_Finalize().
-    void load(const std::string& filepath);
+    // The period column is always validated (finite, > 0); with check_obs the
+    // observed tt and vel columns are too. Aborts all ranks on invalid rows.
+    void load(const std::string& filepath, bool check_obs = false);
+
+    // Log an error listing rows whose value in `col` is non-finite or <= 0.
+    // Returns the number of such rows (0 = column is valid).
+    int check_positive(const real_t* col, const std::string& col_name) const;
 
     // gather synthetic travel times
     void gather_syn_tt();
@@ -134,6 +140,7 @@ class SrcRec {
 
     void get_events();
     void get_periods();
+    std::string filepath_;  // CSV path, for diagnostics
     int nsrc_total_ = 0;
     // Number of rows after load()
     int n_obs_ = 0;