#!/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