Tidied up utils.MinCure and added some specific unit tests. Tried imp…

…roving performance with Numba but it made it slower, so removed.

tests/test_utils.py

@@ -11,11 +11,41 @@
  radius_from_dls,
  get_toolface,
  get_arc,
+ MinCurve
 )
+import json
 
 LAT, LON = (52, 4, 43.1868, "N"), (4, 17, 19.6368, "E")
 
 
+def test_mincurve(decimals=2):
+ with open("tests/test_data/clearance_iscwsa_well_data.json") as f:
+ data = json.load(f)
+
+ for well, survey in data.get('wells').items():
+ def _get_start_xyz(survey):
+ return np.array([
+ survey.get('E')[0],
+ survey.get('N')[0],
+ survey.get('TVD')[0]
+ ])
+
+ mc = MinCurve(
+ md=survey.get('MD'),
+ inc=np.radians(survey.get('IncDeg')),
+ azi=np.radians(survey.get('AziDeg')),
+ start_xyz=_get_start_xyz(survey)
+ )
+
+ assert np.allclose(
+ np.round(mc.poss, decimals), np.round(np.array([
+ survey.get('E'), survey.get('N'), survey.get('TVD')
+ ]).T, decimals)
+ ), "Unexpected position."
+
+ pass
+
+
 def _generate_random_dms(n: int, ndigits: int = None) -> NDArray:
  assert n % 2 == 0, "n must be an even int"
  deg = np.random.randint(0, 180, n)
@@ -235,6 +265,7 @@ def test_get_toolface():
 
 
 def main():
+ test_mincurve()
  test_dms2decimal2dms()
  pass
 

welleng/utils.py

@@ -12,6 +12,10 @@
  NUMBA = False
 
 
+def numbafy(func):
+ func = njit(func)
+
+
 class MinCurve:
  def __init__(
  self,
@@ -39,6 +43,7 @@ def __init__(
  severity.
 
  """
+
  assert unit == "meters" or unit == "feet", (
  'Unknown unit, please select "meters" of "feet"'
  )
@@ -54,98 +59,116 @@ def __init__(
 
  # make two slices with a difference or 1 index to enable array
  # calculations
- md_1 = md[:-1]
- md_2 = md[1:]
-
- inc_1 = inc[:-1]
- inc_2 = inc[1:]
+ md_1, md_2 = self._split_arr(np.array(md))
+ inc_1, inc_2 = self._split_arr(np.array(inc))
+ azi_1, azi_2 = self._split_arr(np.array(azi))
 
- azi_1 = azi[:-1]
- azi_2 = azi[1:]
-
- # calculate the dogleg
- # temp = np.arccos(
- # np.cos(inc_2 - inc_1)
- # - (np.sin(inc_1) * np.sin(inc_2))
- # * (1 - np.cos(azi_2 - azi_1))
- # )
-
- self.dogleg = np.zeros(survey_length)
- self.dogleg[1:] = get_dogleg(inc_1, azi_1, inc_2, azi_2)
+ self.dogleg = self._get_dogleg(survey_length, inc_1, azi_1, inc_2, azi_2)
 
  # calculate rf and assume rf is 1 where dogleg is 0
- self.rf = np.ones(survey_length)
- idx = np.where(self.dogleg != 0)
- self.rf[idx] = 2 / self.dogleg[idx] * np.tan(self.dogleg[idx] / 2)
+ self.rf = self._get_rf(survey_length, self.dogleg)
 
  # calculate the change in md between survey stations
- temp = np.array(md_2) - np.array(md_1)
- self.delta_md = np.zeros(survey_length)
- self.delta_md[1:] = temp
+ self.delta_md = np.diff(self.md, prepend=0)
+
+ args = (
+ self.delta_md, inc_1, azi_1, inc_2, azi_2, self.rf, survey_length
+ )
 
  # calculate change in y direction (north)
- temp = (
- self.delta_md[1:]
+ self.delta_y = self._get_delta_y(*args)
+
+ # calculate change in x direction (east)
+ self.delta_x = self._get_delta_x(*args)
+
+ # calculate change in z direction
+ self.delta_z = self._get_delta_z(*args)
+
+ # calculate the dog leg severity
+ self.dls = self._get_dls(
+ self.dogleg, self.delta_md, survey_length, unit
+ )
+
+ # cumulate the coordinates and add surface coordinates
+ self.poss = np.vstack(
+ np.cumsum(
+ np.array([self.delta_x, self.delta_y, self.delta_z]).T, axis=0
+ ) + self.start_xyz
+ )
+
+ @staticmethod
+ def _get_dogleg(survey_length, inc1, azi1, inc2, azi2):
+ dogleg = np.zeros(survey_length)
+ dogleg[1:] = np.arccos(
+ np.cos(inc2 - inc1)
+ - (np.sin(inc1) * np.sin(inc2))
+ * (1 - np.cos(azi2 - azi1))
+ )
+ return dogleg
+
+ @staticmethod
+ def _split_arr(arr):
+ return (arr[:-1], arr[1:])
+
+ @staticmethod
+ def _get_rf(survey_length, dogleg):
+ rf = np.ones(survey_length)
+ idx = np.where(dogleg != 0)
+ rf[idx] = 2 / dogleg[idx] * np.tan(dogleg[idx] / 2)
+ return rf
+
+ @staticmethod
+ def _get_delta_y(delta_md, inc_1, azi_1, inc_2, azi_2, rf, survey_length):
+ delta_y = np.zeros(survey_length)
+ delta_y[1:] = (
+ delta_md[1:]
  / 2
  * (
  np.sin(inc_1) * np.cos(azi_1)
  + np.sin(inc_2) * np.cos(azi_2)
  )
- * self.rf[1:]
+ * rf[1:]
  )
- self.delta_y = np.zeros(survey_length)
- self.delta_y[1:] = temp
+ return delta_y
 
- # calculate change in x direction (east)
- temp = (
- self.delta_md[1:]
+ @staticmethod
+ def _get_delta_x(delta_md, inc_1, azi_1, inc_2, azi_2, rf, survey_length):
+ delta_x = np.zeros(survey_length)
+ delta_x[1:] = (
+ delta_md[1:]
  / 2
  * (
  np.sin(inc_1) * np.sin(azi_1)
  + np.sin(inc_2) * np.sin(azi_2)
  )
- * self.rf[1:]
+ * rf[1:]
  )
- self.delta_x = np.zeros(survey_length)
- self.delta_x[1:] = temp
+ return delta_x
 
- # calculate change in z direction
- temp = (
- self.delta_md[1:]
+ @staticmethod
+ def _get_delta_z(delta_md, inc_1, azi_1, inc_2, azi_2, rf, survey_length):
+ delta_z = np.zeros(survey_length)
+ delta_z[1:] = (
+ delta_md[1:]
  / 2
  * (np.cos(inc_1) + np.cos(inc_2))
- * self.rf[1:]
+ * rf[1:]
  )
- self.delta_z = np.zeros(survey_length)
- self.delta_z[1:] = temp
-
- # calculate the dog leg severity
+ return delta_z
+
+ @staticmethod
+ def _get_dls(dogleg, delta_md, survey_length, unit):
+ dls = np.zeros(survey_length)
  with np.errstate(divide='ignore', invalid='ignore'):
- temp = np.degrees(self.dogleg[1:]) / self.delta_md[1:]
- self.dls = np.zeros(survey_length)
+ temp = np.degrees(dogleg[1:]) / delta_md[1:]
  mask = np.where(temp != np.nan)
- self.dls[1:][mask] = temp[mask]
+ dls[1:][mask] = temp[mask]
 
  if unit == "meters":
- self.dls *= 30
+ dls *= 30
  else:
- self.dls *= 100
-
- # cumulate the coordinates and add surface coordinates
- self.poss = np.vstack(
- np.cumsum(
- np.array([self.delta_x, self.delta_y, self.delta_z]).T, axis=0
- ) + self.start_xyz
- )
-
-
-def get_dogleg(inc1, azi1, inc2, azi2):
- dogleg = np.arccos(
- np.cos(inc2 - inc1)
- - (np.sin(inc1) * np.sin(inc2))
- * (1 - np.cos(azi2 - azi1))
- )
- return dogleg
+ dls *= 100
+ return dls
 
 
 def get_vec(inc, azi, nev=False, r=1, deg=True):
@@ -228,10 +251,6 @@ def _get_angles(vec):
  return np.stack((inc, azi), axis=1)
 
 
-if NUMBA:
- _get_angles = njit(_get_angles)
-
-
 def get_angles(
  vec: Annotated[NDArray, Literal["N", 3]], nev: bool = False
 ):
@@ -402,7 +421,7 @@ def get_unit_vec(vec):
 
 def linear_convert(data, factor):
  flag = False
- if type(data) != list:
+ if not isinstance(data, list):
  flag = True
  data = [data]
  converted = [d * factor if d is not None else None for d in data]
@@ -513,10 +532,6 @@ def _get_arc_pos_and_vec(dogleg, radius):
  return (pos, vec)
 
 
-if NUMBA:
- _get_arc_pos_and_vec = njit(_get_arc_pos_and_vec)
-
-
 class Arc:
  def __init__(self, dogleg, radius):
  """
@@ -542,18 +557,6 @@ def __init__(self, dogleg, radius):
 
  self.pos, self.vec = _get_arc_pos_and_vec(dogleg, radius)
 
- # self.pos = np.array([
- # np.cos(dogleg),
- # 0.,
- # np.sin(dogleg)
- # ]) * radius
- # self.pos[0] = radius - self.pos[0]
-
- # self.vec = np.array([
- # np.sin(dogleg),
- # 0.,
- # np.cos(dogleg)
- # ])
 
  def transform(self, toolface, pos=None, vec=None, target=False):
  """
@@ -895,3 +898,11 @@ def get_toolface(pos1: NDArray, vec1: NDArray, pos2: NDArray) -> float:
  pos = r.apply(pos2 - pos1)
 
  return np.arctan2(*(np.flip(pos[:2])))
+
+
+if NUMBA:
+ NUMBA_FUNCS = (
+ _get_angles, _get_arc_pos_and_vec
+ )
+ for func in NUMBA_FUNCS:
+ numbafy(func)

welleng/version.py

@@ -1 +1 @@
-__version__ = '0.8.3'
+__version__ = '0.8.4'