All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Pages
frame.h
1 /* Copyright (C) 2004-2014 MBSim Development Team
2  *
3  * This library is free software; you can redistribute it and/or
4  * modify it under the terms of the GNU Lesser General Public
5  * License as published by the Free Software Foundation; either
6  * version 2.1 of the License, or (at your option) any later version.
7  *
8  * This library is distributed in the hope that it will be useful,
9  * but WITHOUT ANY WARRANTY; without even the implied warranty of
10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11  * Lesser General Public License for more details.
12  *
13  * You should have received a copy of the GNU Lesser General Public
14  * License along with this library; if not, write to the Free Software
15  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
16  *
17  * Contact: martin.o.foerg@googlemail.com
18  */
19 
20 #ifndef _FRAME_H__
21 #define _FRAME_H__
22 
23 #include "mbsim/element.h"
24 
25 #ifdef HAVE_OPENMBVCPPINTERFACE
26 #include "mbsim/utils/boost_parameters.h"
27 #include "mbsim/utils/openmbv_utils.h"
28 #endif
29 
30 namespace MBSim {
31 
39  class Frame : public Element {
40  public:
44  enum Feature {
45  position, localPosition, firstTangent, normal, secondTangent, cosy, position_cosy, velocity, angularVelocity, velocity_cosy, velocities, velocities_cosy, angle, dotAngle, all
46  };
47 
52  Frame(const std::string &name = "dummy");
53 
57  virtual ~Frame() {}
58 
59  /* INHERITED INTERFACE ELEMENT */
60  std::string getType() const { return "Frame"; }
61  virtual void plot(double t, double dt = 1);
62  virtual void closePlot();
63  /***************************************************/
64 
65  /* INTERFACE FOR DERIVED CLASSES */
66  virtual int gethSize(int i=0) const { return hSize[i]; }
67  virtual int gethInd(int i=0) const { return hInd[i]; }
68  //virtual ObjectInterface* getParent() { return parent; }
69  //virtual void setParent(ObjectInterface* parent_) { parent = parent_; }
70  virtual const fmatvec::Vec3& getPosition() const { return WrOP; }
71  virtual const fmatvec::Vec3& getLocalPosition() const { return LrOP; }
72  virtual const fmatvec::SqrMat3& getOrientation() const { return AWP; }
73  virtual fmatvec::Vec3& getPosition() { return WrOP; }
74  virtual fmatvec::SqrMat3& getOrientation() { return AWP; }
75  virtual void setPosition(const fmatvec::Vec3 &v) { WrOP = v; }
76  virtual void setLocalPosition(const fmatvec::Vec3 &v) { LrOP = v; }
77  virtual void setOrientation(const fmatvec::SqrMat3 &AWP_) { AWP = AWP_; }
78  virtual const fmatvec::Vec3& getVelocity() const { return WvP; }
79  virtual const fmatvec::Vec3& getAngularVelocity() const { return WomegaP; }
80  virtual const fmatvec::Mat3xV& getJacobianOfTranslation(int j=0) const { return WJP[j]; }
81  virtual const fmatvec::Mat3xV& getJacobianOfRotation(int j=0) const { return WJR[j]; }
82  virtual const fmatvec::Vec3& getGyroscopicAccelerationOfTranslation(int j=0) const { return WjP[j]; }
83  virtual const fmatvec::Vec3& getGyroscopicAccelerationOfRotation(int j=0) const { return WjR[j]; }
84  virtual const fmatvec::Vec3& getAcceleration() const { return WaP; }
85  virtual const fmatvec::Vec3& getAngularAcceleration() const { return WpsiP; }
86  virtual void init(InitStage stage);
87 #ifdef HAVE_OPENMBVCPPINTERFACE
88  BOOST_PARAMETER_MEMBER_FUNCTION( (void), enableOpenMBV, tag, (optional (size,(double),1)(offset,(double),1)(transparency,(double),0))) {
89  OpenMBVFrame ombv(size,offset,"[-1;1;1]",transparency);
90  openMBVFrame=ombv.createOpenMBV();
91  }
92  void setOpenMBVFrame(const boost::shared_ptr<OpenMBV::Frame> &frame) { openMBVFrame = frame; }
93  boost::shared_ptr<OpenMBV::Frame> &getOpenMBVFrame() {return openMBVFrame; }
94 #endif
95  /***************************************************/
96 
97  /* GETTER / SETTER */
98  void sethSize(int size, int i=0) { hSize[i] = size; }
99  void sethInd(int ind, int i=0) { hInd[i] = ind; }
100 
101  fmatvec::Vec3& getVelocity() { return WvP; }
102  fmatvec::Vec3& getAngularVelocity() { return WomegaP; }
103  void setVelocity(const fmatvec::Vec3 &v) { WvP = v; }
104  void setAngularVelocity(const fmatvec::Vec3 &omega) { WomegaP = omega; }
105 
106  void setJacobianOfTranslation(const fmatvec::Mat3xV &WJP_, int j=0) { WJP[j]=WJP_; }
107  void setGyroscopicAccelerationOfTranslation(const fmatvec::Vec3 &WjP_, int j=0) { WjP[j]=WjP_; }
108  void setJacobianOfRotation(const fmatvec::Mat3xV &WJR_, int j=0) { WJR[j]=WJR_; }
109  void setGyroscopicAccelerationOfRotation(const fmatvec::Vec3 &WjR_, int j=0) { WjR[j]=WjR_; }
110  fmatvec::Mat3xV& getJacobianOfTranslation(int j=0) { return WJP[j]; }
111  fmatvec::Mat3xV& getJacobianOfRotation(int j=0) { return WJR[j]; }
112  fmatvec::Vec3& getGyroscopicAccelerationOfTranslation(int j=0) { return WjP[j]; }
113  fmatvec::Vec3& getGyroscopicAccelerationOfRotation(int j=0) { return WjR[j]; }
114  fmatvec::Vec3& getAcceleration() { return WaP; }
115  fmatvec::Vec3& getAngularAcceleration() { return WpsiP; }
116  void setAcceleration(const fmatvec::Vec3 &a) { WaP = a; }
117  void setAngularAcceleration(const fmatvec::Vec3 &psi) { WpsiP = psi; }
118  void setAnglesOfOrientation(const fmatvec::Vec3 &angles_) { angles = angles_; }
119  const fmatvec::Vec3& getAnglesOfOrientation() const { return angles; }
120  void setDotAnglesOfOrientation(const fmatvec::Vec3 &dotAngles_ ) { dotAngles = dotAngles_; }
121  const fmatvec::Vec3& getDotAnglesOfOrientation() const { return dotAngles; }
122 
123  virtual void initializeUsingXML(xercesc::DOMElement *element);
124  virtual xercesc::DOMElement* writeXMLFile(xercesc::DOMNode *element);
125  /***************************************************/
126 
127  protected:
129  // * \brief parent object for plot invocation
130  // */
131  //ObjectInterface* parent;
132 
136  int hSize[2], hInd[2];
137 
141  fmatvec::Vec3 WrOP;
142 
146  fmatvec::Vec3 LrOP;
147 
151  fmatvec::SqrMat3 AWP;
152 
156  fmatvec::Vec3 WvP, WomegaP;
157 
161  fmatvec::Mat3xV WJP[2], WJR[2];
162 
166  fmatvec::Vec3 WjP[2], WjR[2];
167 
171  fmatvec::Vec3 WaP, WpsiP;
172 
176  fmatvec::Vec3 angles, dotAngles;
177 
178 #ifdef HAVE_OPENMBVCPPINTERFACE
179  boost::shared_ptr<OpenMBV::Frame> openMBVFrame;
180 #endif
181  };
182 
187  class FixedRelativeFrame : public Frame {
188 
189  public:
190  FixedRelativeFrame(const std::string &name = "dummy", const fmatvec::Vec3 &r=fmatvec::Vec3(), const fmatvec::SqrMat3 &A=fmatvec::SqrMat3(fmatvec::EYE), const Frame *refFrame=0) : Frame(name), R(refFrame), RrRP(r), ARP(A) {}
191 
192  std::string getType() const { return "FixedRelativeFrame"; }
193 
194  virtual void init(InitStage stage);
195 
196  void setRelativePosition(const fmatvec::Vec3 &r) { RrRP = r; }
197  void setRelativeOrientation(const fmatvec::SqrMat3 &A) { ARP = A; }
198  void setFrameOfReference(const Frame *frame) { R = frame; }
199  void setFrameOfReference(const std::string &frame) { saved_frameOfReference = frame; }
200 
201  const fmatvec::Vec3& getRelativePosition() const { return RrRP; }
202  const fmatvec::SqrMat3& getRelativeOrientation() const { return ARP; }
203  const Frame* getFrameOfReference() const { return R; }
204  const fmatvec::Vec3& getWrRP() const { return WrRP; }
205 
206  void updateRelativePosition() { WrRP = R->getOrientation()*RrRP; }
207  void updatePosition() { updateRelativePosition(); setPosition(R->getPosition() + WrRP); }
208  void updateOrientation() { setOrientation(R->getOrientation()*ARP); }
209  void updateVelocity() { setVelocity(R->getVelocity() + crossProduct(R->getAngularVelocity(), WrRP)); }
210  void updateAngularVelocity() { setAngularVelocity(R->getAngularVelocity()); }
211  void updateStateDependentVariables() {
212  updatePosition();
213  updateOrientation();
214  updateVelocity();
215  updateAngularVelocity();
216  }
217  void updateJacobians(int j=0) {
218  fmatvec::SqrMat3 tWrRP = tilde(WrRP);
219  setJacobianOfTranslation(R->getJacobianOfTranslation(j) - tWrRP*R->getJacobianOfRotation(j),j);
220  setJacobianOfRotation(R->getJacobianOfRotation(j),j);
221  setGyroscopicAccelerationOfTranslation(R->getGyroscopicAccelerationOfTranslation(j) - tWrRP*R->getGyroscopicAccelerationOfRotation(j) + crossProduct(R->getAngularVelocity(),crossProduct(R->getAngularVelocity(),WrRP)),j);
222  setGyroscopicAccelerationOfRotation(R->getGyroscopicAccelerationOfRotation(j),j);
223  }
224  void updateStateDerivativeDependentVariables(const fmatvec::Vec &ud) {
225  setAcceleration(getJacobianOfTranslation()*ud + getGyroscopicAccelerationOfTranslation());
226  setAngularAcceleration(getJacobianOfRotation()*ud + getGyroscopicAccelerationOfRotation());
227  }
228 
229  virtual void initializeUsingXML(xercesc::DOMElement *element);
230  virtual xercesc::DOMElement* writeXMLFile(xercesc::DOMNode *element);
231 
232  protected:
233  const Frame *R;
234  fmatvec::Vec3 RrRP, WrRP;
235  fmatvec::SqrMat3 ARP;
236  std::string saved_frameOfReference;
237  };
238 
239 }
240 
241 #endif /* _FRAME_H_ */
242 
Feature
different interest features for frames
Definition: frame.h:44
fmatvec::Vec3 WaP
acceleration and angular acceleration of coordinate system in inertial frame of reference ...
Definition: frame.h:171
virtual void closePlot()
closes plot file
Definition: frame.cc:90
int hSize[2]
parent object for plot invocation
Definition: frame.h:136
virtual void init(InitStage stage)
plots time series header
Definition: frame.cc:169
cartesian frame on rigid bodies
Definition: frame.h:187
basic class of MBSim mainly for plotting
Definition: element.h:58
std::string getType() const
Definition: frame.h:60
virtual void plot(double t, double dt=1)
plots time dependent data
Definition: frame.cc:48
fmatvec::Vec3 WrOP
position of coordinate system in inertial frame of reference
Definition: frame.h:141
virtual ~Frame()
destructor
Definition: frame.h:57
fmatvec::Mat3xV WJP[2]
Jacobians of translation and rotation from coordinate system to inertial frame.
Definition: frame.h:161
virtual void init(InitStage stage)
plots time series header
Definition: frame.cc:96
InitStage
The stages of the initialization.
Definition: element.h:97
fmatvec::SqrMat3 AWP
transformation matrix in inertial frame of reference
Definition: frame.h:151
std::string name
name of element
Definition: element.h:290
fmatvec::Vec3 WvP
velocity and angular velocity of coordinate system in inertial frame of reference ...
Definition: frame.h:156
cartesian frame on bodies used for application of e.g. links and loads
Definition: frame.h:39
fmatvec::Vec3 LrOP
position of coordinate system in object local frame of reference
Definition: frame.h:146
fmatvec::Vec3 angles
orientation angles and the time derivative of angles of the contour point
Definition: frame.h:176
Vector< Row, AT > crossProduct(const Vector< Row, AT > &x, const Vector< Row, AT > &y)
std::string getType() const
Definition: frame.h:192
SquareMatrix< Row, AT > tilde(const Vector< Row, AT > &x)
fmatvec::Vec3 WjP[2]
Definition: frame.h:166
Frame(const std::string &name="dummy")
constructor
Definition: frame.cc:40

Impressum / Disclaimer / Datenschutz Generated by doxygen 1.8.5 Valid HTML