### Santorini-Amorgos fault-plane workflow MATLAB workflow for identifying, clustering, and visualizing recurrent fault-plane orientations activated during the 2025 Santorini–Amorgos seismic crisis. #### Author and contributions **Code and fault-plane analysis:** Beata Orlecka-Sikora Department of Seismology, Institute of Geophysics, Polish Academy of Sciences (IGF PAS) The fault-plane identification, spatial clustering, PCA plane fitting, orientation clustering, and associated analysis workflow were developed and performed by Beata Orlecka-Sikora. Focal-mechanism calculations were performed separately by Anastasios Kostoglou (University of Thessaloniki) and are not part of the present code package unless explicitly included with the coauthor's agreement. Vasilis Karakostas and Eleftheria Papadimitriou (Aristotle University of Thessaloniki) performed the earthquake relocation and determined the earthquake catalogue parameters, and contributed to the interpretation of the results and manuscript preparation. Code organization, documentation, and repository preparation were assisted by ChatGPT (OpenAI). Scientific design, parameter selection, validation, interpretation, and responsibility for the analysis remain with the author. #### Scope The repository contains the MATLAB code used for the structural analysis of the 2025 Santorini--Amorgos earthquake sequence. The workflow identifies spatial earthquake clusters, fits planar structures to the hypocentres, groups planes according to their orientations, and visualizes the resulting fault-plane families. Two complementary analyses are performed: 1. **Phase-resolved analysis.** Each of the 16 seismic phases is analysed separately. Hypocentres are spatially clustered, planes are fitted to individual clusters, and planes with similar orientations are grouped within each phase. 2. **Complete-catalogue analysis (`PhaseAll`).** The same structural workflow is applied independently to the complete earthquake catalogue without temporal phase segmentation. This analysis yields the six global orientation sets P1--P6 used in the `All Phases` stereonet and Table 1 of the associated manuscript. The phase-resolved and complete-catalogue analyses are therefore independent applications of the same structural workflow. #### Repository files * ##### fault\_plane\_workflow.m Main analysis code. It performs: * conversion of epicentral coordinates to UTM coordinates; * coordinate outlier filtering; * adaptive spatial K-means clustering of hypocentres; * PCA-based plane fitting using singular value decomposition (SVD) to spatial clusters; * evaluation of the RMS distance of hypocentres from each fitted plane; * hierarchical clustering of fitted-plane orientations; * selection of representative planes; * assignment of representative planes to the global P1--P6 orientation families; * export of numerical and graphical results. Strike and dip are reported using the geological right-hand rule (RHR): ``` dip direction = strike + 90 degrees (mod 360) ``` * ##### collect\_dominant\_planes.m Collects the representative planes obtained independently for Phases 1--16 into: ``` PlaneFitStereonet\\\_Dominant\\\_ALL.txt ``` The collector does not recalculate or modify the plane orientations. * ##### plot\_phase\_stereonets.m Produces equal-area stereonet plots of all fitted planes and the representative planes for an individual seismic phase. Example: ``` matlab plot\\\_phase\\\_stereonets('Phase15') ``` The function automatically reads the corresponding plane-fit and grouping files and writes the two stereonet figures. * ##### global\_plane\_clustering.m Performs hierarchical orientation clustering of the representative planes collected from the 16 phase-resolved analyses. This provides a structural consistency check on the six-family architecture obtained from the complete-catalogue analysis. ##### Input data The main workflow requires an earthquake catalogue containing at least: 1. longitude \[degrees]; 2. latitude \[degrees]; 3. depth \[km]; 4. magnitude. Phase-resolved input files correspond to Phases 1--16. The complete catalogue is analysed separately as `PhaseAll`. #### Analysis workflow #### 1\. Spatial coordinate filtering Before K-means clustering, the earthquake coordinates are standardized independently in UTM Easting (`X`), UTM Northing (`Y`), and depth (`Z`) using z-scores. An earthquake is retained only when all three coordinate components satisfy: ``` |zX| <= OUTLIER\\\_Z |zY| <= OUTLIER\\\_Z |zZ| <= OUTLIER\\\_Z ``` For the analyses, the standard threshold is: ``` OUTLIER\\\_Z = 1.5 ``` For Phase 15: ``` OUTLIER\\\_Z = 2.0 ``` `OUTLIER\\\_Z = 1.5` removes an event when at least one of its three spatial coordinates lies more than 1.5 standard deviations from the corresponding coordinate mean. Increasing the threshold to 2.0 makes the filter less restrictive and retains a somewhat broader spatial distribution before clustering and plane fitting. The filter is applied only as an initial spatial screening step; it does not modify the coordinates of retained earthquakes. #### 2\. Spatial K-means clustering K-means is applied to the three-dimensional hypocentral coordinates. The number of spatial clusters is explored adaptively rather than prescribed a priori. The minimum admissible number of earthquakes in a spatial cluster (`MIN\\\_PTS\\\_PER\\\_CLUSTER`) is set to **1, 3, or 5**, depending on the size and spatial distribution of the catalogue for the analysed phase. This parameter prevents the adaptive search from accepting spatial partitions containing clusters smaller than the selected minimum. The phase-resolved analyses use: ``` K-means replicates = 200 ``` The complete-catalogue (`PhaseAll`) analysis uses: ``` K-means replicates = 20 ``` The replicate parameter applies only to K-means clustering of earthquake hypocentres. It is used to reduce dependence on the random initialization of K-means centroids. It does not apply to the subsequent hierarchical clustering of plane orientations. #### 3\. Plane fitting and RMS criterion A plane is fitted to each spatial earthquake cluster using PCA/SVD. The geometric fit is quantified by the root-mean-square distance of hypocentres from the fitted plane. The target geometric coherence is: ``` RMS target = 500 m ``` The 500 m value is a target criterion used in selecting the spatial partition, not an unconditional rejection threshold. If no admissible partition reaches the target within the explored K range, the valid partition minimizing the maximum cluster RMS is retained. #### 4\. Orientation clustering Fitted planes are compared using the three-dimensional angular separation of their normals. Planes are grouped by hierarchical agglomerative clustering using: ``` linkage method = average orientation cutoff = 45 degrees angular penalty = separations > 80 degrees are set to 180 degrees ``` Branches merging below the 45-degree cutoff belong to the same orientation group. #### 5\. Representative planes Within each orientation group, a representative fitted plane is selected as the plane closest to the mean orientation of that group. Representative phase-specific planes are subsequently associated with the global structural families P1--P6. #### 6\. Complete-catalogue analysis The complete earthquake catalogue is analysed independently without division into the 16 temporal phases. The complete-catalogue run uses: ``` K maximum = 30 MIN\\\_PTS\\\_PER\\\_CLUSTER = 5 K-means replicates = 20 OUTLIER\\\_Z = 1.5 RMS target = 500 m orientation cutoff = 45 degrees ``` This analysis independently recovers six dominant orientation sets forming the global P1--P6 structural model. #### Main parameter summary \--- Parameter Phase-resolved analysis Complete catalogue (`PhaseAll`) \--- RMS target 500 m 500 m K maximum 20 30 `MIN\\\_PTS\\\_PER\\\_CLUSTER` 1, 3, or 5 depending on 5 phase K-means replicates 200 20 `OUTLIER\\\_Z` 1.5; Phase 15: 2.0 1.5 RNG seed 1234 1234 Orientation cutoff 45 degrees 45 degrees Angular penalty 80 degrees 80 degrees threshold #### Hierarchical linkage average The small phase-to-phase differences in `MIN\\\_PTS\\\_PER\\\_CLUSTER` and the Phase 15 outlier threshold reflect differences in catalogue size and spatial sampling. They control admissible spatial partitions at the K-means stage; the subsequent plane fitting and orientation-clustering procedure is unchanged. #### Running the analysis For an individual phase, set the corresponding catalogue as the input to `fault\\\_plane\\\_workflow.m` and use the phase-specific settings defined in the code. The workflow generates files containing: * fitted-plane parameters; * spatial cluster assignments; * orientation-group assignments; * representative planes; * run parameters; * diagnostic figures. After all 16 phase-resolved analyses have been completed, run: ``` matlab collect\\\_dominant\\\_planes ``` to collect their representative planes. The complete catalogue is analysed separately with the `PhaseAll` settings listed above. #### Plotting stereonets For example: ``` matlab plot\\\_phase\\\_stereonets('Phase6') ``` The plotting function reads: ``` PlaneFitResults\\\_Phase6.txt StereonetGroups\\\_Phase6.mat ``` and writes: ``` PlaneFitStereonet\\\_AllPlanes\\\_Phase6.png PlaneFitStereonet\\\_DominantPlanes\\\_Phase6.png ``` The stereonet plotting routine preserves the projection convention used for the analysis and figures. #### MATLAB requirements The workflow requires MATLAB and the Statistics and Machine Learning Toolbox. Functions used include: ``` kmeans linkage cluster squareform zscore ``` A latitude/longitude-to-UTM conversion function is also required. If `deg2utm.m` is used, its original authorship and licensing conditions should be retained; it is not automatically covered by the MIT licence applied to the original code in this repository. #### Reproducibility The repository contains the final analysis code and the parameter definitions required to reproduce the workflow. Parameter choices relevant to individual phases are defined directly in the analysis code. #### Citation The citation for the associated manuscript will be added after publication. Until then, users should cite the archived repository version and DOI. #### Licence The original code in this repository is released under the MIT License. See `LICENSE`. Third-party software and external datasets retain their original ownership and licensing conditions.