#!/usr/bin/env python3 # -*- coding: utf-8 -*- """ Created on Thu Aug 11 16:37:16 2022 ################################ Curv2ortho (for LSPIV) ################################ Python library for ortho-rectification of curved surfaces @author: guillaume.bodart """ import numpy as np class polynom(): """ Defines a polynom Contains --------------- - lstCoeff [npArray] Coefficients of the polynom - lim [tuple] """ def __init__(self,coeffs = np.zeros(0), lim = None): if not isinstance(coeffs, np.ndarray): print("\nERROR: coeffs must be a numpy array") return 1 if not (isinstance(lim, tuple) or isinstance(lim, type(None))): print("\nERROR: lim must be a tuple") return 1 self.lstCoeff = coeffs self.lim = lim def set_coeff(self,coeffs,lim): """ Set the coefficient of the polynom Parameters ---------- coeffs : npArray Array containing the coefficients of the polynom sorted in descending order. Returns ------- errno 0 - OK 1 - error with the input. """ if not isinstance(coeffs, np.ndarray): print("\nERROR: coeffs must be a numpy array") return 1 if not isinstance(lim, tuple): print("\nERROR: lim must be a tuple") return 1 # Add coeff self.lstCoeff = coeffs self.lim = lim return 0 def compute_Y(self,X): """ Compute the image of X by the polynomial function Parameters ---------- X : float Entry. Returns ------- out Results. """ if len(self.lstCoeff)!=0: if (isinstance(X,np.ndarray)): # If X is an array or list, then adjust according to the limits if isinstance(self.lim,tuple): X = X[X > self.lim[0]] X = X[X <= self.lim[1]] else: if isinstance(self.lim,tuple): if X < self.lim[0] or X > self.lim[1]: return np.nan out = 0 for ind,el in enumerate(self.lstCoeff[::-1]): out = out + el*X**(ind) return out else: return 0 def compute_Yderive(self,X): """ Compute the image of X by the derivative of the polynomial function Parameters ---------- X : float Entry. Returns ------- out Results. """ if len(self.lstCoeff)!=0: # Compute the derivative lstDeriv = np.zeros(0) for ind,el in enumerate(self.lstCoeff[::-1]): lstDeriv = np.append(lstDeriv,el*ind) lstDeriv = lstDeriv[1:] out = 0 for ind,el in enumerate(lstDeriv): out = out + el*X**(ind) return out else: return 0 def lenAB(poly,X): """ Compute the length of the polynomial function 'f' within the x sequence. The function is defined as y=f(x). The step corresponds to the x sampling distance Parameters ---------- poly : polynome() polynomial function describing the curve of the flow surface. X : npArray x sequence used. Returns ------- float length of the function whithin the given x sequence. """ s = 0 step = np.mean(X[1:]-X[:-1]) for el in X: s = s + np.sqrt(1+poly.compute_Yderive(el)**2)*step return np.abs(s) def computeLongMesh(poly,X,nbV): """ Compute the longitudinal mesh according to the polynomial equation given and the x sequence. nbV corresponds to the number of points along the curve (i.e. the number of pixels on the i axis of the ortho-rectified image) Parameters ---------- poly : polynome polynomial function describing the curve of the flow surface. X : npArray x sequence considered. nbV : int Number of points along the curve (i.e. the number of pixels on the i axis of the ortho-rectified image) . Returns ------- xy : TYPE DESCRIPTION. """ # Create position vector xy = np.empty((nbV,2)) # Compute the mesh step size d = lenAB(poly,X)/nbV # Compute first point xy[0,:] = X[0],poly.compute_Y(X[0]) k = 1 for el in xy[:-1]: # Compute vector director xD = 1 yD = poly.compute_Yderive(el[0]) nD = np.sqrt(xD**2+yD**2) xD = xD/nD yD = yD/nD # Compute position x2 x2 = el[0]+(d*xD) y2 = el[1]+(d*yD) xy[k,:] = x2,y2 k = k+1 return xy class function(): """ Defines a function made of polynoms Contains --------------- - poly [polynom()] List of polynom() (Class containing the polynom) - lim [tuple] List containing the limits of each polynom """ def __init__(self): self.poly = [] self.lim = [] def add_poly(self,poly): if not isinstance(poly,polynom): print("\nERROR: poly must be a polynom()\n") return 1 self.poly.append(poly) self.lim.append(poly.lim) return 0 def compute_Y(self,X): if isinstance(X,np.ndarray): out = np.empty(0) for el in self.poly: out = np.concatenate((out,el.compute_Y(X))) return out else: for el in self.poly: if np.isfinite(el.compute_Y(X)): out = el.compute_Y(X) return out def lenAB(self,X): """ Compute the length of the polynomial function 'f' within the x sequence. The function is defined as y=f(x). The step corresponds to the x sampling distance Parameters ---------- poly : polynome() polynomial function describing the curve of the flow surface. X : npArray x sequence used. Returns ------- float length of the function whithin the given x sequence. """ s = 0 step = np.mean(X[1:]-X[:-1]) for poly in self.poly: for el in X: if np.isfinite(poly.compute_Yderive(el)): #check out nan values s = s + np.sqrt(1+poly.compute_Yderive(el)**2)*step return np.abs(s) def computeLongMesh(self,X,nbV): """ Compute the longitudinal mesh according to the polynomial equation given and the x sequence. nbV corresponds to the number of points along the curve (i.e. the number of pixels on the i axis of the ortho-rectified image) Parameters ---------- poly : polynome polynomial function describing the curve of the flow surface. X : npArray x sequence considered. nbV : int Number of points along the curve (i.e. the number of pixels on the i axis of the ortho-rectified image) . Returns ------- xy : TYPE DESCRIPTION. """ # Create position vector xy = np.empty((nbV,2)) # Compute the mesh step size d = self.lenAB(X)/nbV # Compute first point xy[0,:] = X[0],self.compute_Y(X[0]) k = 1 for el in xy[:-1]: # Compute vector director xD = 1 yD = self.compute_Yderive(el[0]) nD = np.sqrt(xD**2+yD**2) xD = xD/nD yD = yD/nD # Compute position x2 x2 = el[0]+(d*xD) y2 = el[1]+(d*yD) xy[k,:] = x2,y2 k = k+1 return xy