Merge pull request #23 from TomoATT/devel

v0.2.4

.github/workflows/build-docs.yml

@@ -2,10 +2,11 @@
 # For more information see: https://docs.github.com/en/actions/automating-builds-and-tests/building-and-testing-python
 
 name: Build documentations
-
+permissions:
+  contents: write
 on:
   push:
-    branches: [ "docs" ]
+    branches: [ "docs", "main" ]
   workflow_run:
     workflows: ["Upload Python Package"]
     types:

.github/workflows/build-test-conda.yml

@@ -1,8 +1,11 @@
 name: Python Package using Conda
-
+permissions:
+  contents: write
 on:
   push:
     branches: [ devel ]
+  pull_request:
+    branches: [ devel, main ]
 
 jobs:
   build-linux:

.github/workflows/python-publish.yml

@@ -7,14 +7,12 @@
 # documentation.
 
 name: Upload Python Package
-
+permissions:
+  contents: write
 on:
-  pull_request:
+  push:
     branches: [ main ]
 
-permissions:
-  contents: read
-
 jobs:
   deploy:
 
@@ -37,11 +35,5 @@ jobs:
       with:
         user: __token__
         password: ${{ secrets.PYPI_API_TOKEN }}
-
-  # Build_docs:
-  #   needs: deploy
-  #   uses: MIGG-NTU/PyTomoATT/.github/workflows/build-docs.yml@devel
-  #   with:
-  #     config-path: .github/workflows/build-docs.yml
-  #   secrets: inherit
+
 

pyproject.toml

@@ -6,7 +6,7 @@ build-backend = "hatchling.build"
 name = "PyTomoATT"
 description = "Your project description"
 authors = [
-    { name = "Mijian Xu", email = "[email protected]" },
+    { name = "Mijian Xu", email = "[email protected]" },
     { name = "Masaru Nagaso" }
 ]
 license = { text = "GPL-3.0" }

pytomoatt/_version.py

@@ -1 +1 @@
-__version__ = '0.2.3'
+__version__ = '0.2.4'

pytomoatt/para.py

@@ -2,6 +2,7 @@
 from .utils.common import init_axis, str2val
 
 yaml = YAML()
+yaml.default_flow_style = True
 
 class ATTPara:
     """Class for read and write parameter file with ``yaml`` format
@@ -24,10 +25,10 @@ def init_axis(self):
             self.input_params['domain']['n_rtp'],
         )
         return dep, lat, lon, dd, dt, dp
-    
+
     def update_param(self, key: str, value) -> None:
         """Update a parameter in the YAML file.
-        
+
         :param key: The key of parameter file to be set. Use '.' to separate the keys.
         :type key: str
         """

pytomoatt/src_rec.py

@@ -3,7 +3,7 @@
 import pandas as pd
 from .distaz import DistAZ
 from .setuplog import SetupLog
-from .utils.src_rec_utils import define_rec_cols, setup_rec_points_dd, get_rec_points_types
+from .utils.src_rec_utils import define_rec_cols, setup_rec_points_dd, get_rec_points_types, update_position
 from sklearn.metrics.pairwise import haversine_distances
 import copy
 
@@ -220,8 +220,8 @@ def read(cls, fname: str, dist_in_data=False, name_net_and_sta=False, **kwargs):
         """
         sr = cls(fname=fname, **kwargs)
         alldf = pd.read_table(
-            fname, sep="\s+", header=None, comment="#"
-        )
+                fname, sep=r"\s+", header=None, comment="#", low_memory=False
+            )
 
         last_col_src = 12
         dd_col = 11
@@ -1489,7 +1489,7 @@ def add_noise(self, range_in_sec=0.1, mean_in_sec=0.0, shape="gaussian"):
                 )
             rec_type["tt"] += noise
 
-    def rotate(self, clat:float, clon:float, angle:float):
+    def rotate(self, clat:float, clon:float, angle:float, reverse=False):
         """Rotate sources and receivers around a center point
 
         :param clat: Latitude of the center
@@ -1499,36 +1499,19 @@ def rotate(self, clat:float, clon:float, angle:float):
         :param angle: anti-clockwise angle in degree
         :type angle: float
         """
-        from .utils.rotate import rtp_rotation
+        from .utils.rotate import rtp_rotation, rtp_rotation_reverse
+
+        rotation_func = rtp_rotation_reverse if reverse else rtp_rotation
 
-        self.sources["evla"], self.sources["evlo"] = rtp_rotation(
-            self.sources["evla"], self.sources["evlo"], clat, clon, angle
+        self.sources["evla"], self.sources["evlo"] = rotation_func(
+            self.sources["evla"].to_numpy(), self.sources["evlo"].to_numpy(), clat, clon, angle
         )
-        for i, row in self.src_points.iterrows():
-            self.src_points.loc[i, "evla"], self.src_points.loc[i, "evlo"] = rtp_rotation(
-                row["evla"], row["evlo"], clat, clon, angle
-            )
 
-        self.receivers["stla"], self.receivers["stlo"] = rtp_rotation(
-            self.receivers["stla"], self.receivers["stlo"], clat, clon, angle
+        self.receivers["stla"], self.receivers["stlo"] = rotation_func(
+            self.receivers["stla"].to_numpy(), self.receivers["stlo"].to_numpy(), clat, clon, angle
         )
-        for i, row in self.rec_points.iterrows():
-            self.rec_points.loc[i, "stla"] = self.receivers[self.receivers["staname"] == row["staname"]]["stla"]
-            self.rec_points.loc[i, "stlo"] = self.receivers[self.receivers["staname"] == row["staname"]]["stlo"]
-
-        if not self.rec_points_cs.empty:
-            for i, row in self.rec_points_cs.iterrows():
-                self.rec_points_cs.loc[i, "stla1"] = self.receivers[self.receivers["staname"] == row["staname1"]]["stla"]
-                self.rec_points_cs.loc[i, "stlo1"] = self.receivers[self.receivers["staname"] == row["staname1"]]["stlo"]
-                self.rec_points_cs.loc[i, "stla2"] = self.receivers[self.receivers["staname"] == row["staname2"]]["stla"]
-                self.rec_points_cs.loc[i, "stlo2"] = self.receivers[self.receivers["staname"] == row["staname2"]]["stlo"]
-
-        if not self.rec_points_cr.empty:
-            for i, row in self.rec_points_cr.iterrows():
-                self.rec_points_cr.loc[i, "stla"] = self.receivers[self.receivers["staname"] == row["staname"]]["stla"]
-                self.rec_points_cr.loc[i, "stlo"] = self.receivers[self.receivers["staname"] == row["staname"]]["stlo"]
-                self.rec_points_cr.loc[i, "evla2"] = self.sources[self.sources["event_id"] == row["event_id2"]]["evla"]
-                self.rec_points_cr.loc[i, "evlo2"] = self.sources[self.sources["event_id"] == row["event_id2"]]["evlo"]
+
+        update_position(self)
 
     def to_utm(self, zone):
         """Convert sources and receivers to UTM coordinates
@@ -1547,68 +1530,8 @@ def to_utm(self, zone):
             self.receivers["stlo"], self.receivers["stla"]
         )
 
-        self.src_points = self.src_points.merge(
-            self.sources[['event_id', 'evlo', 'evla']],
-            on='event_id',
-            how='left',
-            suffixes=('', '_new')
-        )
-        self.src_points['evlo'] = self.src_points['evlo_new']
-        self.src_points['evla'] = self.src_points['evla_new']
-        self.src_points.drop(columns=['evlo_new', 'evla_new'], inplace=True)
-
-        self.rec_points = self.rec_points.merge(
-            self.receivers[['staname', 'stlo', 'stla']],
-            on='staname',
-            how='left',
-            suffixes=('', '_new')
-        )
-        self.rec_points['stlo'] = self.rec_points['stlo_new']
-        self.rec_points['stla'] = self.rec_points['stla_new']
-        self.rec_points.drop(columns=['stlo_new', 'stla_new'], inplace=True)
-
-        if not self.rec_points_cs.empty:
-            self.rec_points_cs = self.rec_points_cs.merge(
-                self.receivers[['staname', 'stlo', 'stla']],
-                left_on='staname1',
-                right_on='staname',
-                how='left',
-            )
-            self.rec_points_cs['stlo1'] = self.rec_points_cs['stlo']
-            self.rec_points_cs['stla1'] = self.rec_points_cs['stla']
-            self.rec_points_cs.drop(columns=['stlo', 'stla', 'staname'], inplace=True)
-
-            self.rec_points_cs = self.rec_points_cs.merge(
-                self.receivers[['staname', 'stlo', 'stla']],
-                left_on='staname2',
-                right_on='staname',
-                how='left',
-            )
-            self.rec_points_cs['stlo2'] = self.rec_points_cs['stlo']
-            self.rec_points_cs['stla2'] = self.rec_points_cs['stla']
-            self.rec_points_cs.drop(columns=['stlo', 'stla', 'staname'], inplace=True)
-
-        if not self.rec_points_cr.empty:
-            self.rec_points_cr = self.rec_points_cr.merge(
-                self.receivers[['staname', 'stlo', 'stla']],
-                on='staname',
-                how='left',
-                suffixes=('', '_new')
-            )
-            self.rec_points_cr['stlo'] = self.rec_points_cr['stlo_new']
-            self.rec_points_cr['stla'] = self.rec_points_cr['stla_new']
-            self.rec_points_cr.drop(columns=['stlo_new', 'stla_new'], inplace=True)
-
-            self.rec_points_cr = self.rec_points_cr.merge(
-                self.sources[['event_id', 'evlo', 'evla']],
-                left_on='event_id2',
-                right_on='event_id',
-                how='left',
-            )
-            self.rec_points_cr['evlo2'] = self.rec_points_cr['evlo']
-            self.rec_points_cr['evla2'] = self.rec_points_cr['evla']
-            self.rec_points_cr.drop(columns=['evlo', 'evla', 'event_id'], inplace=True)
-
+        update_position(self)
+
     def write_receivers(self, fname: str):
         """
         Write receivers to a txt file.

pytomoatt/utils/common.py

@@ -136,21 +136,38 @@ def ignore_nan_3d(data):
     xidx = np.arange(data.shape[2])
     zz, xx, yy = np.meshgrid(zidx, yidx, xidx, indexing='ij')
     interpolated = griddata(
-        points, values, 
-        (zz, xx, yy), 
+        points, values,
+        (zz, xx, yy),
         method='nearest'
     )
     result = interpolated.reshape(data.shape)
     return result
 
 def str2val(str_val):
+    # single value handling
+    # return integer
     try:
         return int(str_val)
     except ValueError:
-        try:
-            return float(str_val)
-        except ValueError:
-            try:
-                return [float(v) for v in str_val.strip('[]').split(',')]
-            except ValueError:
-                return str_val
+        pass
+
+    # return float
+    try:
+        return float(str_val)
+    except ValueError:
+        pass
+
+    # list values handling
+    # return list of integer
+    try:
+        return [int(v) for v in str_val.strip('[]').split(',')]
+    except ValueError:
+        pass
+
+    # return list of float
+    try:
+        return [float(v) for v in str_val.strip('[]').split(',')]
+    except ValueError:
+        pass
+
+    return str_val

pytomoatt/utils/rotate.py

@@ -11,39 +11,51 @@ def rtp2xyz(r,theta,phi):
     return (x,y,z)
 
 # Cartesian coordinates to Spherical coordinate
-def xyz2rtp(x,y,z):
-    # theta: -90~90;  phi: -180~180
-    r = np.sqrt(x**2+y**2+z**2)
-    theta = np.arcsin(z/r)
-    phi = np.arcsin(y/r/np.cos(theta))
-    idx = np.where((phi > 0) & (x*y < 0))
-    if len(idx[0]) > 0:
-        phi[idx] = np.pi - phi[idx]
-    idx = np.where((phi < 0) & (x*y > 0))
-    if len(idx[0]) > 0:
-        phi[idx] = -np.pi - phi[idx]
-    return (r,theta*RAD2DEG,phi*RAD2DEG)
+def xyz2rtp(x, y, z):
+    """Convert Cartesian coordinates (x, y, z) to spherical coordinates (r, theta, phi).
+
+    Args:
+        x (float or np.ndarray): X coordinate(s).
+        y (float or np.ndarray): Y coordinate(s).
+        z (float or np.ndarray): Z coordinate(s).
+
+    Returns:
+        tuple: A tuple containing radius (r), polar angle (theta), and azimuthal angle (phi) in degrees.
+    """
+    r = np.sqrt(x**2 + y**2 + z**2)
+
+    # theta = arctan(z / sqrt(x^2 + y^2))
+    theta = np.arctan2(z, np.sqrt(x**2 + y**2))
+
+    # phi = arctan(y / x) 
+    phi = np.arctan2(y, x)
+
+    # Convert radians to degrees
+    theta = theta * RAD2DEG
+    phi = phi * RAD2DEG
+
+    return r, theta, phi
 
 # anti-clockwise rotation along x-axis
 def rotate_x(x,y,z,theta):
     new_x = x 
     new_y = y *  np.cos(theta*DEG2RAD) + z * -np.sin(theta*DEG2RAD)
     new_z = y *  np.sin(theta*DEG2RAD) + z *  np.cos(theta*DEG2RAD)
-    return (new_x,new_y,new_z)
+    return new_x, new_y, new_z
 
 # anti-clockwise rotation along y-axis
 def rotate_y(x,y,z,theta):
     new_x = x *  np.cos(theta*DEG2RAD) + z *  np.sin(theta*DEG2RAD)
     new_y = y
     new_z = x * -np.sin(theta*DEG2RAD) + z *  np.cos(theta*DEG2RAD)
-    return (new_x,new_y,new_z)
+    return new_x, new_y, new_z
 
 # anti-clockwise rotation along z-axis
 def rotate_z(x,y,z,theta):
     new_x = x *  np.cos(theta*DEG2RAD) + y * -np.sin(theta*DEG2RAD)
     new_y = x *  np.sin(theta*DEG2RAD) + y *  np.cos(theta*DEG2RAD)
     new_z = z 
-    return (new_x,new_y,new_z)
+    return new_x, new_y, new_z
 
 
 # rotate to the new coordinate, satisfying the center r0,t0,p0 -> r0,0,0 and a anticlockwise angle psi
@@ -61,9 +73,9 @@ def rtp_rotation(t,p,theta0,phi0,psi):
     (x,y,z) = rotate_x(x,y,z,psi)
 
     # step 5: x,y,z -> r,t,p
-    (new_r,new_t,new_p) = xyz2rtp(x,y,z)
+    _, new_t, new_p = xyz2rtp(x,y,z)
 
-    return (new_t,new_p)
+    return new_t, new_p
 
 
 def rtp_rotation_reverse(new_t,new_p,theta0,phi0,psi):
@@ -80,6 +92,6 @@ def rtp_rotation_reverse(new_t,new_p,theta0,phi0,psi):
     (x,y,z) = rotate_z(x,y,z,phi0)
 
     # step 5: x,y,z -> r,t,p
-    (r,t,p) = xyz2rtp(x,y,z)
+    _, t, p = xyz2rtp(x,y,z)
 
-    return (t,p)
+    return t, p