timeseries: output velocities in units per year (#341)

* `timeseries`: output velocities in units per year

* [skip ci] Add note on per-year scaling

src/dolphin/timeseries.py

@@ -375,7 +375,7 @@ def estimate_velocity_pixel(x: ArrayLike, y: ArrayLike, w: ArrayLike) -> Array:
     Returns
     -------
     velocity : np.array, 0D
-        The estimated velocity in (unw unit) / day.
+        The estimated velocity in (unw unit) / year.
 
     """
     # Jax polyfit will grab the first *2* dimensions of y to solve in a batch
@@ -402,8 +402,8 @@ def estimate_velocity(
     Returns
     -------
     velocity : np.array 2D
-        The estimated velocity in (unw unit) / day.
-        E.g. if the unwrapped phase is in radians, the velocity is in rad/day.
+        The estimated velocity in (unw unit) / year calculated as 365.25 * rad/day.
+        E.g. if the unwrapped phase is in radians, the velocity is in rad/year.
 
     """
     # TODO: weighted least squares using correlation?
@@ -427,7 +427,9 @@ def estimate_velocity(
         velos = vmap(estimate_velocity_pixel, in_axes=(None, -1, -1))(
             x_arr, unw_pixels, weights_pixels
         )
-    return velos.reshape(n_rows, n_cols)
+    # Currently `velos` is in units / day,
+    days_per_year = 365.25
+    return velos.reshape(n_rows, n_cols) * days_per_year
 
 
 def datetime_to_float(dates: Sequence[DateOrDatetime]) -> np.ndarray:
@@ -466,6 +468,10 @@ def create_velocity(
 ) -> None:
     """Perform pixel-wise (weighted) linear regression to estimate velocity.
 
+    The units of `output_file` are in (unwrapped units) / year.
+    E.g. if the files in `unw_file_list` are in radians, the output velocity
+    is in radians / year, which is calculated as 365.25 * radians / day.
+
     Parameters
     ----------
     unw_file_list : Sequence[PathOrStr]

tests/test_timeseries.py

@@ -13,7 +13,8 @@
 DT = 12
 START_DATE = datetime(2020, 1, 1)
 SHAPE = 50, 50
-VELO_RAD_PER_DAY = 0.2  # rad / day
+VELO_RAD_PER_YEAR = 5.0
+VELO_RAD_PER_DAY = VELO_RAD_PER_YEAR / 365.25
 
 
 def make_sar_dates():
@@ -203,7 +204,7 @@ def expected_velo(self, data, x_arr):
         for i in range(SHAPE[0]):
             for j in range(SHAPE[1]):
                 out[i, j] = np.polyfit(x_arr, sar_phases[:, i, j], 1)[0]
-        return out
+        return out * 365.25
 
     def test_stack(self, data, x_arr, expected_velo):
         sar_dates, sar_phases, ifg_date_pairs, ifgs = data