HOW TO RUN AND TEST THE CODE

WHAT YOU NEED (one-time setup)
1. Install Python 3 from https://www.python.org/downloads/ (tick "Add Python to PATH"
   on Windows during install).
2. Open a terminal:
     Windows: open "Command Prompt"
     Mac:     open "Terminal"
3. Install the required libraries by typing this one line and pressing Enter:
     pip install numpy matplotlib scipy pandas

THE TWO GROUPS OF SCRIPTS
There are two groups. The first runs out of the box and reproduces the synthetic
comparison model that validates the method. The second reproduces the headline
result from the real Ryugu boulder catalogue and needs one extra data file.

GROUP A - SYNTHETIC MODEL (runs immediately, no extra data)
4. Put all the code files in one folder, go to that folder in the terminal:
     Windows example:  cd C:\Users\YourName\Downloads\Code
     Mac example:      cd ~/Downloads/Code
5. Run each script with "python" then the file name:
     python code_monte_carlo.py
     python code_hessian_stability.py
     python code_deff_sensitivity.py
     python code_trap_number.py
     python code_generate_figures.py

   Expected results (synthetic model):
     code_monte_carlo.py       Ryugu about 0.44 +/- 0.05, Bennu about 0.038,
                               Itokawa about 0.051
     code_hessian_stability.py two eigenvalues, about -9.5e-05 and +1.0e-05,
                               labelled SADDLE POINT
     code_deff_sensitivity.py  T_max of about 1.15, 0.42, 0.19, 0.11 at
                               30, 50, 75, 100 m  (SYNTHETIC model trend)
     code_trap_number.py       a table of T_max for the three asteroids
     code_generate_figures.py  creates Fig1.png to Fig6.png in the same folder

GROUP B - REAL RYUGU CATALOGUE (reproduces the paper's headline numbers)
6. The headline Ryugu result uses the real boulder catalogue of Ray et al. (2026),
   which is those authors' dataset and is cited, not included here. Download it from
   the figshare deposit listed in Boulder_Data_Sources.md, save it as
   Ryugu_Boulder_Attributes.csv in the same folder, then run:
     python code_real_ryugu_trap_number.py

   Expected results (real catalogue):
     T_max about 0.5 at Otohime (0.34 shape-corrected to 0.67 long-axis),
     and the d_eff sensitivity 1.9, 0.67, 0.31, 0.19 at 30, 50, 75, 100 m.
   These 1.9 / 0.67 / 0.31 / 0.19 values are the ones quoted in the manuscript;
   the synthetic code_deff_sensitivity.py trend in Group A is the model analogue.
     python run_trap_number_real_catalogue.py  is a generic version that accepts
     any longitude/latitude/size CSV.

NOTES
- The analysis scripts only print numbers to the screen; they create no files
  (except code_generate_figures.py, which writes the figure images).
- Synthetic numbers vary slightly each run because the boulder field is random,
  except where a fixed starting seed is set; the averages match the values above.
- The Group A analysis scripts read the core file code_trap_number.py, so keep all
  files together in one folder.

------------------------------------------------------------------
REAL RYUGU CATALOGUE (required by code_real_ryugu_trap_number.py and
verify_from_scratch.py)
  File: Ryugu_Boulder_Attributes.csv (included in this Code/ folder
  and in ../Data/). Source: Ray, A., Ruj, T., Komatsu, G., Vincent,
  L.M., Kikuchi, H., Hemmi, R., Ormo, J. (2026), Scientific Reports
  16, 14404, https://doi.org/10.1038/s41598-026-40370-4 (CC BY 4.0).
  Columns: FID, Longitude_X, Latitude_Y, Length_m, Orientation, Aspect.

INDEPENDENT VERIFICATION (added)
  verify_from_scratch.py     - recomputes g_surface, T_max (0.67 long
                               axis / 0.34 shape-corrected), the d_eff
                               sweep, the Eq 5 = Eq 6 identity, and a
                               dimensional check, all from the real
                               catalogue, printing PASS/FAIL per check.
  verify_hessian_and_plot.py - Hessian stability (saddle) at Otohime
                               and regenerates the trap-number map.
  Both run with: python <script>.py   (needs numpy, scipy, pandas,
  matplotlib). They reproduce the manuscript values from first
  principles using the simple form T = G*M_local/(g_surface*d_eff^2),
  i.e. the steepest-slope limit sin(alpha)=1.
