## 1 Introduction  
The advance of Global Navigation Satellite System (GNSS) technology is a valuable asset for the study of crustal deformation, providing insights to the strain and velocity fields and, furthermore, to the geodynamic evolution of a region (e.g., [1, 2]). Crustal deformation is the result of elastic strain accumulation and, when a threshold is exceeded, almost all of the elastic deformation turns permanent by fault rupturing [3] and, consequently, into earthquake genesis. Mapping and measuring strain rates by combining GNSS and Interferometric SAR provides not only a significant temporal resolution, but, most importantly, a high spatial resolution and coverage, which are essential to fault-specific seismic hazard assessment [4]. The strain rate can be calculated by applying geological (e.g., [5]), seismological (e.g., [6]), and geodetic methods (e.g., [7, 8]), or a combination of these [9]. During the last few decades, the latter has mostly been based on the continuous monitoring of a permanent GNSS station network and/or Interferometric Synthetic Aperture Radar (InSAR), providing time series [10] which are well suited to study interseismic strain accumulation (e.g., [11]).  
Satellite geodetic analysis is a useful tool for estimating the deformation of the upper crust in a regional scale. Although it covers a relatively short-term period of the seismic cycle (a few years up to a few decades) that maybe inadequate for capturing the full length (which can range from a few hundred up to several thousands of years), it provides a valuable extended snapshot of the present-day deformation field. As a result, it must be emphasized that the estimated short-term geodetic rate may not necessarily reflect the long-term deformation rate, as the latter may not be stable during the seismic cycle [12; 13]. Moreover, no conclusion can be drawn about the maturity stage of the fault, and, hence, the time of the next rupture. Nevertheless, satellite geodesy can depict the strain accumulation rate along known faults or along faults which were never mapped before with conventional techniques or provide important strain rate thresholds.  
Our study area consists of the Strymon basin in North Greece, which develops in a NW-SE direction from the Greek-Bulgarian border down to the estuary in the Orphanos Gulf (North Aegean Sea). The basin is bounded by active faults controlling the marginal topography, causing subsidence, while low strain rates (30-70 \(\times\) 10\({}^{-9}\)/yr) are documented by broad-coverage geodetic studies [7; 14; 15; 16; 17]. Low strain rates infer long seismic recurrence intervals; as a result, many faults have not been activated during the historically and instrumentally recorded periods in earthquake catalogues. Moreover, the extracted geological rates are elusive or subtle due to low strain rate per fault, and seismological data are insufficient due to the long recurrence intervals and the incompleteness of the historical earthquake record. Hence, geodetic data are valuable in such a tectonic regime.  
The geodetic datasets' processing is based on the triangulation methodology (its basic principle is based on combinations of three different GNSS stations), while the extracted results are further geost statistically processed, implementing interpolation (kriging) methodology, in order to determine the results' spatial distribution in a grid pattern.  
Taking all the above into account, the scope of this study is to: (i) estimate the upper crust deformation, based on primary geodetic data, (ii) offer a comparison between GNNS and InSAR results, and (iii) correlate geodetically extracted strain with the surrounding active faults and define the potential seismic hazard.