Commit 98e36730 authored by Michele Maris's avatar Michele Maris
Browse files

u

parent 524eb845
Loading
Loading
Loading
Loading

CHANGELOG

deleted100644 → 0
+0 −6
Original line number Diff line number Diff line
2021 dec 01 - Repo created
2021 dec 02 - Uploaded PyPRM 0.0 
2021 dec 03 - Added example
2023 jul 12 - PyPRM 1.1.0 : Updated wobble angles operator in the original it was rotZ(omegaVA)rotY(zVAX) in the new is rotZ(omegaVA)rotY(zVAX)rotZ(-omegaVA)
2023 jul 14 - V1.1.1 : Update structure to allow creation of a PyPI package
2023 jul 18 - V1.1.2 : Updated version to allow comments

CHANGELOG.md

0 → 100644
+20 −0
Original line number Diff line number Diff line
# Changelog

<!--next-version-placeholder-->

## v0.0.0 (01/12/2021)

- Repo created!

## v0.0.1 (02/12/2021)

- Uploaded PyPRM 0.0!

## v0.0.2 (03/12/2021)

- Added example!

# 2023 jul 12 - PyPRM 1.1.0 : Updated wobble angles operator in the original it was rotZ(omegaVA)rotY(zVAX) in the new is rotZ(omegaVA)rotY(zVAX)rotZ(-omegaVA)
# 2023 jul 14 - V1.1.1 : Update structure to allow creation of a PyPI package
# 2023 jul 18 - V1.1.2 : Updated version to allow comments
+8 −4
Original line number Diff line number Diff line
@@ -10,11 +10,11 @@ $ pip install lspe_prm

## Usage

'''PYTHON
```python
import numpy as np
from lspe_prm import PRM0

# instantiate object as an ideal telescope setup
## instantiate object as an ideal telescope setup
prm=PRM0(ideal=True)

# set some angles
@@ -33,9 +33,13 @@ theta_phi=np.array([
#creates the list of pointing matrices
RotMatr=prm(theta_phi)

#to recover pointings makes dot product with (0,0,1) vector
#to compute pointings use the dot product with (0,0,1) vector
Pointing=np.array([k.dot(np.array([0,0,1])) for k in RotMatr]).T
'''

#to compute orientation use the dot product with (1,0,0) vector
Orient=np.array([k.dot(np.array([1,0,0])) for k in RotMatr]).T

```

## Contributing

+1 −366
Original line number Diff line number Diff line
@@ -9,369 +9,4 @@ except :

from .roto_translation import *


__VERSION__="1.1.1 - 2021 Dec 03 - 2023 jul 14"
__DESCRIPTION__="""
M.Maris

order 0 Pointing Reconstruction Model (PRM)

=========================
#Example:

import numpy as np
from lspe_prm import PRM0

# instantiate object as an ideal telescope setup
prm=PRM0(ideal=True)

# set some angles
prm.tiltFor = 0.1 # deg
prm.cphi0   = -5.1 # deg
prm.commit() # this is mandatory

# creates a list of [theta, phi] couples (deg)
theta_phi=np.array([
 [15., 20]
 ,[25., 50]
 ,[35., 30]
 ,[10., 210]
])

#creates the list of pointing matrices
RotMatr=prm(theta_phi)

#to recover pointings makes dot product with (0,0,1) vector
Pointing=np.array([k.dot(np.array([0,0,1])) for k in RotMatr]).T
=========================

Elemental Rotation Matrices according to Aereonautical convention

Those matrices projects Coordinate Axis of the Moving Reference Frame (example: camera) into the fixed reference frame (example: geographic topocentric coordinates).

alpha > 0 for an anticlockwise rotation

anticlockwise means:

Z rotation: 
     Y^
      |
   Z (.)->X 

Z exits from (X,Y) plane, X rotated toward Y 
equivalent to Z(X->Y)

Y rotation: Y exits from (Z,X) plane, Z rotated toward X
     X^
      |
   Y (.)->Z 
equivalent to Y(Z->X)

X rotation: X exits from (Y,Z) plane, Y rotated toward Z
     Z^
      |
   X (.)->Y 
equivalent to X(Y->Z)

RotX(alpha) : 
   1  0  0
   0  c -s
   0  s  c

RotY(alpha) : 
   c  0  s
   0  1  0
  -s  0  c

RotZ(alpha) : 
   c  s  0  
  -s  c  0
   0  0  1

with c=cos(alpha), s=sin(alpha)

Issue 
0.0.0 - 2021 Oct 23 : First implementation
0.0.1 - 2021 Nov 20 : Definition of matrices in agreement with STRIP PRM Doc
0.0.2 - 2021 Nov 22 : Improved definition of matrices
0.0.3 - 2021 Dec  2 : Structure improved, added example
1.0.0 - 2021 Dec  3 : First Issue after review with the Pointing Working Group
1.1.0 - 2023 Jul 12 : Changed the definition of wobble angles operator
1.1.1 - 2023 Jul 12 : package structure adapted to PyPI requirements
"""

import numpy as np

class PRM0 :
   @property
   def roll(self) :
      return self._d['roll']
   @roll.setter
   def roll(self,this) :
      self._d['roll']=this
      self._mark['roll']=1
      self._matr['roll']=self.ERotZ(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def pan(self) :
      return self._d['pan']
   @pan.setter
   def pan(self,this) :
      self._d['pan']=this
      self._mark['pan']=1
      self._matr['pan']=self.ERotY(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def tilt(self) :
      return self._d['tilt']
   @tilt.setter
   def tilt(self,this) :
      self._d['tilt']=this
      self._mark['tilt']=1
      self._matr['tilt']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def ctheta0(self) :
      return self._d['ctheta0']
   @ctheta0.setter
   def ctheta0(self,this) :
      self._d['ctheta0']=this
      self._mark['ctheta0']=1
      self._matr['ctheta0']=self.ERotY(-np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def cphi0(self) :
      return self._d['cphi0']
   @cphi0.setter
   def cphi0(self,this) :
      self._d['cphi0']=this
      self._mark['cphi0']=1
      self._matr['cphi0']=self.ERotZ(-np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def tiltFork(self) :
      return self._d['tiltFork']
   @tiltFork.setter
   def tiltFork(self,this) :
      self._d['tiltFork']=this
      self._mark['tiltFork']=1
      self._matr['tiltFork']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def zVAX(self) :
      return self._d['zVAX']
   @zVAX.setter
   def zVAX(self,this) :
      self._d['zVAX']=this
      self._mark['zVAX']=1
      self._matr['zVAX']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def omegaVAX(self) :
      return self._d['omegaVAX']
   @omegaVAX.setter
   def omegaVAX(self,this) :
      self._d['omegaVAX']=this
      self._mark['omegaVAX']=1
      self._matr['omegaVAX']=self.ERotZ(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def isComplete(self) :
      return self._isComplete
   #
   @property
   def isAllSet(self) :
      for k in self.keys() :
         if self._mark[k]==0 :
            return False
      return True
   #
   def __init__(self,ideal=False) :
      self._isComplete=False
      self._d={}
      for k in self.keys() :
         self._d[k]=None
      self._mark={}
      for k in self.keys() :
         self._mark[k]=0
      self._matr={}
      for k in self.keys() :
         self._matr[k]=None
      if ideal :
         self.setIdeal()
   #
   def keys(self) :
      return ['roll','pan','tilt','ctheta0','cphi0','tiltFork','zVAX','omegaVAX']
   #
   def __getitem__(self,this) :
      return self._d[this]
   #
   def __setitem__(self,this,that) :
      self._d[this]=that
   #
   def M(self,k) :
      return self._matr[k]
   #
   def copy(self) :
      import copy
      return copy.deepcopy(self)
   #
   def _base_matrix(self,Angle) :
      _angle=np.array([Angle]) if np.isscalar(Angle) else Angle
      c=np.cos(_angle)
      s=np.sin(_angle)
      I=np.zeros([len(_angle),3,3])
      I[:,0,0]=1
      I[:,1,1]=1
      I[:,2,2]=1
      return I,c,s
   #
   def ERotX(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,1,1]= c
      R[:,1,2]=-s
      R[:,2,1]= s
      R[:,2,2]= c
      return R
   #
   def ERotY(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,0,0]= c
      R[:,0,2]= s
      R[:,2,0]=-s
      R[:,2,2]= c
      return R
   #
   def ERotZ(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,0,0]= c
      R[:,0,1]=-s
      R[:,1,0]= s
      R[:,1,1]= c
      return R
   #
   def setIdeal(self) :
      self.roll=0.
      self.pan=0.
      self.tilt=0.
      self.tiltFork=0.
      self.zVAX=0.
      self.omegaVAX=0.
      self.cphi0=0.
      self.ctheta0=0.
      self.commit()
   #
   def commit(self,oldWobble=False) :
      """ if oldWobble==True the old version of wobble angle operator is applied 
      
          current version : M3 = matr['omegaVAX'].dot(matr['zVAX'].dot(matr['omegaVAX'].T))
          old version     : M3 = matr['omegaVAX'].dot(matr['zVAX'])
          
          the old version does not handle properly the cases zVAX=0.
      """
      if self.isAllSet :
         out=self._matr['tilt'][0].dot(self._matr['pan'][0].dot(self._matr['roll'][0]))
         out=self._matr['ctheta0'][0].dot(out)
         out.shape=(1,3,3)
         self._matr['M1']=out*1
         #
         out=self._matr['cphi0'][0].dot(self._matr['tiltFork'][0])
         out.shape=(1,3,3)
         self._matr['M2']=out*1
         #
         if oldWobble :
            #old version
            out=self._matr['omegaVAX'][0].dot(self._matr['zVAX'][0])
         else :
            #new version
            out=self._matr['omegaVAX'][0].dot(self._matr['zVAX'][0].dot(self._matr['omegaVAX'][0].T))
         out.shape=(1,3,3)
         self._matr['M3']=out*1
         #
         self._isComplete=True
      else :
         self._isComplete=False
         raise Error('Configuration not completed','')
   #
   def __call__(self,CTHETA_PHI) :
      """ 
Apply rotations for a list of NPOS control angles: (CTR_THETA, CTR_PHI)

CTHETA_PHI.shape = (NPOS,2)
          
CTHETA_PHI is an array of NPOS positions, each position being defined by the couple of angles (CTR_THETA,CTR_PHI) expressed in deg.
          
The output is a (3,3,NPOS) array, representing the corresponding NPOS matrices.
          
Example:
          
CTHETA_PHI=np.array([
   [45.,0]
  ,[45.,90]
  ,[45.,180]
  ,[45.,270]
])

PMatr=prm(CTHETA_PHI)

print(PMatr.shape) 

will result in (3,3,4)

while PMatr[:,:,0] is the rotation matrix for the case (45.,0.) and PMatr[:,2,0] is the corresponding pointing vector.

A different ordering could be to have (NPOS,3,3) matrices, this can be obtained by using the np.moveaxis function

PMatrReorder=np.moveaxis(PMatr,2,-1)

Now the first matrix is PMatrReorder[0] and the corresponding pointing vector is PMatrReorder[0][2].

the option axis0isNPOS=True could be used 

In this case 

PMatrB=prm(CTHETA_PHI,axis0isNPOS=True)

print(PMatrB.shape) 

will result in (4,3,3)

while PMatrB[0] is the rotation matrix for the case (45.,0.) and PMatrB[:,2] is the corresponding pointing vector.

      """ 
      if not self.isComplete :
         raise Error('Configuration not completed or not committed','')
      #
      self._d['CTHETA_PHI']=CTHETA_PHI
      #
      # theta rotation
      self._matr['CTRL_THETA']=self.ERotY(np.deg2rad(CTHETA_PHI[:,0]))
      #
      # phi rotation
      self._matr['CTRL_PHI']=self.ERotZ(np.deg2rad(CTHETA_PHI[:,1]))
      #
      # apply each theta rotation to M1
      out=[k.dot(self._matr['M1'][0]) for k in self._matr['CTRL_THETA']]
      #
      # apply M2 to each out
      out1=[self._matr['M2'][0].dot(k) for k in out]
      #
      # apply phi rotation to each out2
      out2=[self._matr['CTRL_PHI'][ik].dot(k) for ik, k in enumerate(out1)]
      #
      # apply M3 to each out3
      out3=np.array([self._matr['M3'][0].dot(k) for k in out2])
      #
      #stores Rout 
      self._matr['Rout']=out3
      #
      return out3
      
from .lspe_prm import *
+372 −0
Original line number Diff line number Diff line
# read version from installed package



from .roto_translation import *


__VERSION__="1.1.1 - 2021 Dec 03 - 2023 jul 14"
__DESCRIPTION__="""
M.Maris

order 0 Pointing Reconstruction Model (PRM)

=========================
#Example:

import numpy as np
from lspe_prm import PRM0

# instantiate object as an ideal telescope setup
prm=PRM0(ideal=True)

# set some angles
prm.tiltFor = 0.1 # deg
prm.cphi0   = -5.1 # deg
prm.commit() # this is mandatory

# creates a list of [theta, phi] couples (deg)
theta_phi=np.array([
 [15., 20]
 ,[25., 50]
 ,[35., 30]
 ,[10., 210]
])

#creates the list of pointing matrices
RotMatr=prm(theta_phi)

#to recover pointings makes dot product with (0,0,1) vector
Pointing=np.array([k.dot(np.array([0,0,1])) for k in RotMatr]).T
=========================

Elemental Rotation Matrices according to Aereonautical convention

Those matrices projects Coordinate Axis of the Moving Reference Frame (example: camera) into the fixed reference frame (example: geographic topocentric coordinates).

alpha > 0 for an anticlockwise rotation

anticlockwise means:

Z rotation: 
     Y^
      |
   Z (.)->X 

Z exits from (X,Y) plane, X rotated toward Y 
equivalent to Z(X->Y)

Y rotation: Y exits from (Z,X) plane, Z rotated toward X
     X^
      |
   Y (.)->Z 
equivalent to Y(Z->X)

X rotation: X exits from (Y,Z) plane, Y rotated toward Z
     Z^
      |
   X (.)->Y 
equivalent to X(Y->Z)

RotX(alpha) : 
   1  0  0
   0  c -s
   0  s  c

RotY(alpha) : 
   c  0  s
   0  1  0
  -s  0  c

RotZ(alpha) : 
   c  s  0  
  -s  c  0
   0  0  1

with c=cos(alpha), s=sin(alpha)

Issue 
0.0.0 - 2021 Oct 23 : First implementation
0.0.1 - 2021 Nov 20 : Definition of matrices in agreement with STRIP PRM Doc
0.0.2 - 2021 Nov 22 : Improved definition of matrices
0.0.3 - 2021 Dec  2 : Structure improved, added example
1.0.0 - 2021 Dec  3 : First Issue after review with the Pointing Working Group
1.1.0 - 2023 Jul 12 : Changed the definition of wobble angles operator
1.1.1 - 2023 Jul 12 : package structure adapted to PyPI requirements
"""

import numpy as np

class PRM0 :
   @property
   def roll(self) :
      return self._d['roll']
   @roll.setter
   def roll(self,this) :
      self._d['roll']=this
      self._mark['roll']=1
      self._matr['roll']=self.ERotZ(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def pan(self) :
      return self._d['pan']
   @pan.setter
   def pan(self,this) :
      self._d['pan']=this
      self._mark['pan']=1
      self._matr['pan']=self.ERotY(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def tilt(self) :
      return self._d['tilt']
   @tilt.setter
   def tilt(self,this) :
      self._d['tilt']=this
      self._mark['tilt']=1
      self._matr['tilt']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def ctheta0(self) :
      return self._d['ctheta0']
   @ctheta0.setter
   def ctheta0(self,this) :
      self._d['ctheta0']=this
      self._mark['ctheta0']=1
      self._matr['ctheta0']=self.ERotY(-np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def cphi0(self) :
      return self._d['cphi0']
   @cphi0.setter
   def cphi0(self,this) :
      self._d['cphi0']=this
      self._mark['cphi0']=1
      self._matr['cphi0']=self.ERotZ(-np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def tiltFork(self) :
      return self._d['tiltFork']
   @tiltFork.setter
   def tiltFork(self,this) :
      self._d['tiltFork']=this
      self._mark['tiltFork']=1
      self._matr['tiltFork']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def zVAX(self) :
      return self._d['zVAX']
   @zVAX.setter
   def zVAX(self,this) :
      self._d['zVAX']=this
      self._mark['zVAX']=1
      self._matr['zVAX']=self.ERotX(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def omegaVAX(self) :
      return self._d['omegaVAX']
   @omegaVAX.setter
   def omegaVAX(self,this) :
      self._d['omegaVAX']=this
      self._mark['omegaVAX']=1
      self._matr['omegaVAX']=self.ERotZ(np.deg2rad(np.array([this])))
      self._isComplete=False
   #
   @property
   def isComplete(self) :
      return self._isComplete
   #
   @property
   def isAllSet(self) :
      for k in self.keys() :
         if self._mark[k]==0 :
            return False
      return True
   #
   def __init__(self,ideal=False) :
      self._isComplete=False
      self._d={}
      for k in self.keys() :
         self._d[k]=None
      self._mark={}
      for k in self.keys() :
         self._mark[k]=0
      self._matr={}
      for k in self.keys() :
         self._matr[k]=None
      if ideal :
         self.setIdeal()
   #
   def keys(self) :
      return ['roll','pan','tilt','ctheta0','cphi0','tiltFork','zVAX','omegaVAX']
   #
   def __getitem__(self,this) :
      return self._d[this]
   #
   def __setitem__(self,this,that) :
      self._d[this]=that
   #
   def M(self,k) :
      return self._matr[k]
   #
   def copy(self) :
      import copy
      return copy.deepcopy(self)
   #
   def _base_matrix(self,Angle) :
      _angle=np.array([Angle]) if np.isscalar(Angle) else Angle
      c=np.cos(_angle)
      s=np.sin(_angle)
      I=np.zeros([len(_angle),3,3])
      I[:,0,0]=1
      I[:,1,1]=1
      I[:,2,2]=1
      return I,c,s
   #
   def ERotX(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,1,1]= c
      R[:,1,2]=-s
      R[:,2,1]= s
      R[:,2,2]= c
      return R
   #
   def ERotY(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,0,0]= c
      R[:,0,2]= s
      R[:,2,0]=-s
      R[:,2,2]= c
      return R
   #
   def ERotZ(self,Angle) :
      R,c,s=self._base_matrix(Angle)
      R[:,0,0]= c
      R[:,0,1]=-s
      R[:,1,0]= s
      R[:,1,1]= c
      return R
   #
   def setIdeal(self) :
      self.roll=0.
      self.pan=0.
      self.tilt=0.
      self.tiltFork=0.
      self.zVAX=0.
      self.omegaVAX=0.
      self.cphi0=0.
      self.ctheta0=0.
      self.commit()
   #
   def commit(self,oldWobble=False) :
      """ if oldWobble==True the old version of wobble angle operator is applied 
      
          current version : M3 = matr['omegaVAX'].dot(matr['zVAX'].dot(matr['omegaVAX'].T))
          old version     : M3 = matr['omegaVAX'].dot(matr['zVAX'])
          
          the old version does not handle properly the cases zVAX=0.
      """
      if self.isAllSet :
         out=self._matr['tilt'][0].dot(self._matr['pan'][0].dot(self._matr['roll'][0]))
         out=self._matr['ctheta0'][0].dot(out)
         out.shape=(1,3,3)
         self._matr['M1']=out*1
         #
         out=self._matr['cphi0'][0].dot(self._matr['tiltFork'][0])
         out.shape=(1,3,3)
         self._matr['M2']=out*1
         #
         if oldWobble :
            #old version
            out=self._matr['omegaVAX'][0].dot(self._matr['zVAX'][0])
         else :
            #new version
            out=self._matr['omegaVAX'][0].dot(self._matr['zVAX'][0].dot(self._matr['omegaVAX'][0].T))
         out.shape=(1,3,3)
         self._matr['M3']=out*1
         #
         self._isComplete=True
      else :
         self._isComplete=False
         raise Error('Configuration not completed','')
   #
   def __call__(self,CTHETA_PHI) :
      """ 
Apply rotations for a list of NPOS control angles: (CTR_THETA, CTR_PHI)

CTHETA_PHI.shape = (NPOS,2)
          
CTHETA_PHI is an array of NPOS positions, each position being defined by the couple of angles (CTR_THETA,CTR_PHI) expressed in deg.
          
The output is a (3,3,NPOS) array, representing the corresponding NPOS matrices.
          
Example:
          
CTHETA_PHI=np.array([
   [45.,0]
  ,[45.,90]
  ,[45.,180]
  ,[45.,270]
])

PMatr=prm(CTHETA_PHI)

print(PMatr.shape) 

will result in (3,3,4)

while PMatr[:,:,0] is the rotation matrix for the case (45.,0.) and PMatr[:,2,0] is the corresponding pointing vector.

A different ordering could be to have (NPOS,3,3) matrices, this can be obtained by using the np.moveaxis function

PMatrReorder=np.moveaxis(PMatr,2,-1)

Now the first matrix is PMatrReorder[0] and the corresponding pointing vector is PMatrReorder[0][2].

the option axis0isNPOS=True could be used 

In this case 

PMatrB=prm(CTHETA_PHI,axis0isNPOS=True)

print(PMatrB.shape) 

will result in (4,3,3)

while PMatrB[0] is the rotation matrix for the case (45.,0.) and PMatrB[:,2] is the corresponding pointing vector.

      """ 
      if not self.isComplete :
         raise Error('Configuration not completed or not committed','')
      #
      self._d['CTHETA_PHI']=CTHETA_PHI
      #
      # theta rotation
      self._matr['CTRL_THETA']=self.ERotY(np.deg2rad(CTHETA_PHI[:,0]))
      #
      # phi rotation
      self._matr['CTRL_PHI']=self.ERotZ(np.deg2rad(CTHETA_PHI[:,1]))
      #
      # apply each theta rotation to M1
      out=[k.dot(self._matr['M1'][0]) for k in self._matr['CTRL_THETA']]
      #
      # apply M2 to each out
      out1=[self._matr['M2'][0].dot(k) for k in out]
      #
      # apply phi rotation to each out2
      out2=[self._matr['CTRL_PHI'][ik].dot(k) for ik, k in enumerate(out1)]
      #
      # apply M3 to each out3
      out3=np.array([self._matr['M3'][0].dot(k) for k in out2])
      #
      #stores Rout 
      self._matr['Rout']=out3
      #
      return out3