PageRenderTime 36ms CodeModel.GetById 1ms RepoModel.GetById 0ms app.codeStats 0ms

/src/ckw05.c

https://github.com/mattbornski/spice
C | 1112 lines | 260 code | 280 blank | 572 comment | 48 complexity | 30550d8534755a98b2b9b19300e7c33a MD5 | raw file
  1. /* ckw05.f -- translated by f2c (version 19980913).
  2. You must link the resulting object file with the libraries:
  3. -lf2c -lm (in that order)
  4. */
  5. #include "f2c.h"
  6. /* Table of constant values */
  7. static integer c__4 = 4;
  8. static integer c__15 = 15;
  9. static integer c__2 = 2;
  10. static integer c__6 = 6;
  11. static integer c__1 = 1;
  12. /* $Procedure CKW05 ( Write CK segment, type 5 ) */
  13. /* Subroutine */ int ckw05_(integer *handle, integer *subtyp, integer *degree,
  14. doublereal *begtim, doublereal *endtim, integer *inst, char *ref,
  15. logical *avflag, char *segid, integer *n, doublereal *sclkdp,
  16. doublereal *packts, doublereal *rate, integer *nints, doublereal *
  17. starts, ftnlen ref_len, ftnlen segid_len)
  18. {
  19. /* System generated locals */
  20. integer i__1, i__2;
  21. doublereal d__1;
  22. /* Local variables */
  23. integer addr__, i__;
  24. extern /* Subroutine */ int chkin_(char *, ftnlen), dafps_(integer *,
  25. integer *, doublereal *, integer *, doublereal *);
  26. doublereal descr[5];
  27. extern /* Subroutine */ int errch_(char *, char *, ftnlen, ftnlen),
  28. errdp_(char *, doublereal *, ftnlen), dafada_(doublereal *,
  29. integer *);
  30. doublereal dc[2];
  31. extern /* Subroutine */ int dafbna_(integer *, doublereal *, char *,
  32. ftnlen);
  33. integer ic[6];
  34. extern /* Subroutine */ int dafena_(void);
  35. extern logical failed_(void);
  36. integer chrcod, refcod;
  37. extern integer bsrchd_(doublereal *, integer *, doublereal *);
  38. extern /* Subroutine */ int namfrm_(char *, integer *, ftnlen);
  39. extern integer lastnb_(char *, ftnlen);
  40. integer packsz;
  41. extern /* Subroutine */ int sigerr_(char *, ftnlen), chkout_(char *,
  42. ftnlen), setmsg_(char *, ftnlen), errint_(char *, integer *,
  43. ftnlen);
  44. extern integer lstltd_(doublereal *, integer *, doublereal *);
  45. extern logical vzerog_(doublereal *, integer *), return_(void);
  46. integer winsiz;
  47. extern logical odd_(integer *);
  48. /* $ Abstract */
  49. /* Write a type 5 segment to a CK file. */
  50. /* $ Disclaimer */
  51. /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
  52. /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
  53. /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
  54. /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
  55. /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
  56. /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
  57. /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
  58. /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
  59. /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
  60. /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
  61. /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
  62. /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
  63. /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
  64. /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
  65. /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
  66. /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
  67. /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
  68. /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
  69. /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
  70. /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
  71. /* $ Required_Reading */
  72. /* CK */
  73. /* NAIF_IDS */
  74. /* ROTATION */
  75. /* TIME */
  76. /* $ Keywords */
  77. /* POINTING */
  78. /* FILES */
  79. /* $ Declarations */
  80. /* $ Abstract */
  81. /* Declare parameters specific to CK type 05. */
  82. /* $ Disclaimer */
  83. /* THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE */
  84. /* CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S. */
  85. /* GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE */
  86. /* ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE */
  87. /* PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS" */
  88. /* TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY */
  89. /* WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A */
  90. /* PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC */
  91. /* SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE */
  92. /* SOFTWARE AND RELATED MATERIALS, HOWEVER USED. */
  93. /* IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA */
  94. /* BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT */
  95. /* LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND, */
  96. /* INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS, */
  97. /* REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE */
  98. /* REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY. */
  99. /* RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF */
  100. /* THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY */
  101. /* CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE */
  102. /* ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE. */
  103. /* $ Required_Reading */
  104. /* CK */
  105. /* $ Keywords */
  106. /* CK */
  107. /* $ Restrictions */
  108. /* None. */
  109. /* $ Author_and_Institution */
  110. /* N.J. Bachman (JPL) */
  111. /* $ Literature_References */
  112. /* None. */
  113. /* $ Version */
  114. /* - SPICELIB Version 1.0.0, 20-AUG-2002 (NJB) */
  115. /* -& */
  116. /* CK type 5 subtype codes: */
  117. /* Subtype 0: Hermite interpolation, 8-element packets. Quaternion */
  118. /* and quaternion derivatives only, no angular velocity */
  119. /* vector provided. Quaternion elements are listed */
  120. /* first, followed by derivatives. Angular velocity is */
  121. /* derived from the quaternions and quaternion */
  122. /* derivatives. */
  123. /* Subtype 1: Lagrange interpolation, 4-element packets. Quaternion */
  124. /* only. Angular velocity is derived by differentiating */
  125. /* the interpolating polynomials. */
  126. /* Subtype 2: Hermite interpolation, 14-element packets. */
  127. /* Quaternion and angular angular velocity vector, as */
  128. /* well as derivatives of each, are provided. The */
  129. /* quaternion comes first, then quaternion derivatives, */
  130. /* then angular velocity and its derivatives. */
  131. /* Subtype 3: Lagrange interpolation, 7-element packets. Quaternion */
  132. /* and angular velocity vector provided. The quaternion */
  133. /* comes first. */
  134. /* Packet sizes associated with the various subtypes: */
  135. /* End of file ck05.inc. */
  136. /* $ Brief_I/O */
  137. /* Variable I/O Description */
  138. /* -------- --- -------------------------------------------------- */
  139. /* HANDLE I Handle of an CK file open for writing. */
  140. /* SUBTYP I CK type 5 subtype code. */
  141. /* DEGREE I Degree of interpolating polynomials. */
  142. /* BEGTIM I Start time of interval covered by segment. */
  143. /* ENDTIM I End time of interval covered by segment. */
  144. /* INST I NAIF code for a s/c instrument or structure. */
  145. /* REF I Reference frame name. */
  146. /* AVFLAG I True if the segment will contain angular velocity. */
  147. /* SEGID I Segment identifier. */
  148. /* N I Number of packets. */
  149. /* SCLKDP I Encoded SCLK times. */
  150. /* PACKTS I Array of packets. */
  151. /* RATE I Nominal SCLK rate in seconds per tick. */
  152. /* NINTS I Number of intervals. */
  153. /* STARTS I Encoded SCLK interval start times. */
  154. /* MAXDEG P Maximum allowed degree of interpolating polynomial. */
  155. /* $ Detailed_Input */
  156. /* HANDLE is the file handle of a CK file that has been */
  157. /* opened for writing. */
  158. /* SUBTYP is an integer code indicating the subtype of the */
  159. /* the segment to be created. */
  160. /* DEGREE is the degree of the polynomials used to */
  161. /* interpolate the quaternions contained in the input */
  162. /* packets. All components of the quaternions are */
  163. /* interpolated by polynomials of fixed degree. */
  164. /* BEGTIM, */
  165. /* ENDTIM are the beginning and ending encoded SCLK times */
  166. /* for which the segment provides pointing */
  167. /* information. BEGTIM must be less than or equal to */
  168. /* ENDTIM, and at least one data packet must have a */
  169. /* time tag T such that */
  170. /* BEGTIM < T < ENDTIM */
  171. /* - - */
  172. /* INST is the NAIF integer code for the instrument or */
  173. /* structure for which a segment is to be created. */
  174. /* REF is the NAIF name for a reference frame relative to */
  175. /* which the pointing information for INST is */
  176. /* specified. */
  177. /* AVFLAG is a logical flag which indicates whether or not */
  178. /* the segment will contain angular velocity. */
  179. /* SEGID is the segment identifier. A CK segment */
  180. /* identifier may contain up to 40 characters. */
  181. /* N is the number of packets in the input packet */
  182. /* array. */
  183. /* SCLKDP are the encoded spacecraft clock times associated */
  184. /* with each pointing instance. These times must be */
  185. /* strictly increasing. */
  186. /* PACKTS contains a time-ordered array of data packets */
  187. /* representing the orientation of INST relative to */
  188. /* the frame REF. Each packet contains a SPICE-style */
  189. /* quaternion and optionally, depending on the */
  190. /* segment subtype, attitude derivative data, from */
  191. /* which a C-matrix and an angular velocity vector */
  192. /* may be derived. */
  193. /* See the discussion of quaternion styles in */
  194. /* Particulars below. */
  195. /* The C-matrix represented by the Ith data packet is */
  196. /* a rotation matrix that transforms the components */
  197. /* of a vector expressed in the base frame specified */
  198. /* by REF to components expressed in the instrument */
  199. /* fixed frame at the time SCLKDP(I). */
  200. /* Thus, if a vector V has components x, y, z in the */
  201. /* base frame, then V has components x', y', z' */
  202. /* in the instrument fixed frame where: */
  203. /* [ x' ] [ ] [ x ] */
  204. /* | y' | = | CMAT | | y | */
  205. /* [ z' ] [ ] [ z ] */
  206. /* The attitude derivative information in PACKTS(I) */
  207. /* gives the angular velocity of the instrument fixed */
  208. /* frame at time SCLKDP(I) with respect to the */
  209. /* reference frame specified by REF. */
  210. /* The direction of an angular velocity vector gives */
  211. /* the right-handed axis about which the instrument */
  212. /* fixed reference frame is rotating. The magnitude */
  213. /* of the vector is the magnitude of the */
  214. /* instantaneous velocity of the rotation, in radians */
  215. /* per second. */
  216. /* Packet contents and the corresponding */
  217. /* interpolation methods depend on the segment */
  218. /* subtype, and are as follows: */
  219. /* Subtype 0: Hermite interpolation, 8-element */
  220. /* packets. Quaternion and quaternion */
  221. /* derivatives only, no angular */
  222. /* velocity vector provided. */
  223. /* Quaternion elements are listed */
  224. /* first, followed by derivatives. */
  225. /* Angular velocity is derived from */
  226. /* the quaternions and quaternion */
  227. /* derivatives. */
  228. /* Subtype 1: Lagrange interpolation, 4-element */
  229. /* packets. Quaternion only. Angular */
  230. /* velocity is derived by */
  231. /* differentiating the interpolating */
  232. /* polynomials. */
  233. /* Subtype 2: Hermite interpolation, 14-element */
  234. /* packets. Quaternion and angular */
  235. /* angular velocity vector, as well as */
  236. /* derivatives of each, are provided. */
  237. /* The quaternion comes first, then */
  238. /* quaternion derivatives, then */
  239. /* angular velocity and its */
  240. /* derivatives. */
  241. /* Subtype 3: Lagrange interpolation, 7-element */
  242. /* packets. Quaternion and angular */
  243. /* velocity vector provided. The */
  244. /* quaternion comes first. */
  245. /* Angular velocity is always specified relative to */
  246. /* the base frame. */
  247. /* RATE is the nominal rate of the spacecraft clock */
  248. /* associated with INST. Units are seconds per */
  249. /* tick. RATE is used to scale angular velocity */
  250. /* to radians/second. */
  251. /* NINTS is the number of intervals that the pointing */
  252. /* instances are partitioned into. */
  253. /* STARTS are the start times of each of the interpolation */
  254. /* intervals. These times must be strictly increasing */
  255. /* and must coincide with times for which the segment */
  256. /* contains pointing. */
  257. /* $ Detailed_Output */
  258. /* None. See $Particulars for a description of the effect of this */
  259. /* routine. */
  260. /* $ Parameters */
  261. /* MAXDEG is the maximum allowed degree of the interpolating */
  262. /* polynomial. If the value of MAXDEG is increased, */
  263. /* the SPICELIB routine CKPFS must be changed */
  264. /* accordingly. In particular, the size of the */
  265. /* record passed to CKRnn and CKEnn must be */
  266. /* increased, and comments describing the record size */
  267. /* must be changed. */
  268. /* $ Exceptions */
  269. /* If any of the following exceptions occur, this routine will return */
  270. /* without creating a new segment. */
  271. /* 1) If HANDLE is not the handle of a C-kernel opened for writing */
  272. /* the error will be diagnosed by routines called by this */
  273. /* routine. */
  274. /* 2) If the last non-blank character of SEGID occurs past index 40, */
  275. /* the error SPICE(SEGIDTOOLONG) is signaled. */
  276. /* 3) If SEGID contains any nonprintable characters, the error */
  277. /* SPICE(NONPRINTABLECHARS) is signaled. */
  278. /* 4) If the first encoded SCLK time is negative then the error */
  279. /* SPICE(INVALIDSCLKTIME) is signaled. If any subsequent times */
  280. /* are negative the error will be detected in exception (5). */
  281. /* 5) If the encoded SCLK times are not strictly increasing, */
  282. /* the error SPICE(TIMESOUTOFORDER) is signaled. */
  283. /* 6) If the name of the reference frame is not one of those */
  284. /* supported by the routine FRAMEX, the error */
  285. /* SPICE(INVALIDREFFRAME) is signaled. */
  286. /* 7) If the number of packets N is not at least 1, the error */
  287. /* SPICE(TOOFEWPACKETS) will be signaled. */
  288. /* 8) If NINTS, the number of interpolation intervals, is less than */
  289. /* or equal to 0, the error SPICE(INVALIDNUMINTS) is signaled. */
  290. /* 9) If the encoded SCLK interval start times are not strictly */
  291. /* increasing, the error SPICE(TIMESOUTOFORDER) is signaled. */
  292. /* 10) If an interval start time does not coincide with a time for */
  293. /* which there is an actual pointing instance in the segment, */
  294. /* then the error SPICE(INVALIDSTARTTIME) is signaled. */
  295. /* 11) This routine assumes that the rotation between adjacent */
  296. /* quaternions that are stored in the same interval has a */
  297. /* rotation angle of THETA radians, where */
  298. /* 0 < THETA < pi. */
  299. /* _ */
  300. /* The routines that evaluate the data in the segment produced */
  301. /* by this routine cannot distinguish between rotations of THETA */
  302. /* radians, where THETA is in the interval [0, pi), and */
  303. /* rotations of */
  304. /* THETA + 2 * k * pi */
  305. /* radians, where k is any integer. These "large" rotations will */
  306. /* yield invalid results when interpolated. You must ensure that */
  307. /* the data stored in the segment will not be subject to this */
  308. /* sort of ambiguity. */
  309. /* 12) If any quaternion has magnitude zero, the error */
  310. /* SPICE(ZEROQUATERNION) is signaled. */
  311. /* 13) If the interpolation window size implied by DEGREE is not */
  312. /* even, the error SPICE(INVALIDDEGREE) is signaled. The window */
  313. /* size is DEGREE+1 for Lagrange subtypes and is (DEGREE+1)/2 */
  314. /* for Hermite subtypes. */
  315. /* 14) If an unrecognized subtype code is supplied, the error */
  316. /* SPICE(NOTSUPPORTED) is signaled. */
  317. /* 15) If DEGREE is not at least 1 or is greater than MAXDEG, the */
  318. /* error SPICE(INVALIDDEGREE) is signaled. */
  319. /* 16) If the segment descriptor bounds are out of order, the */
  320. /* error SPICE(BADDESCRTIMES) is signaled. */
  321. /* 17) If there is no element of SCLKDP that lies between BEGTIM and */
  322. /* ENDTIM inclusive, the error SPICE(EMPTYSEGMENT) is signaled. */
  323. /* 18) If RATE is zero, the error SPICE(INVALIDVALUE) is signaled. */
  324. /* $ Files */
  325. /* A new type 5 CK segment is written to the CK file attached */
  326. /* to HANDLE. */
  327. /* $ Particulars */
  328. /* This routine writes a CK type 5 data segment to the open CK */
  329. /* file according to the format described in the type 5 section of */
  330. /* the CK Required Reading. The CK file must have been opened with */
  331. /* write access. */
  332. /* Quaternion Styles */
  333. /* ----------------- */
  334. /* There are different "styles" of quaternions used in */
  335. /* science and engineering applications. Quaternion styles */
  336. /* are characterized by */
  337. /* - The order of quaternion elements */
  338. /* - The quaternion multiplication formula */
  339. /* - The convention for associating quaternions */
  340. /* with rotation matrices */
  341. /* Two of the commonly used styles are */
  342. /* - "SPICE" */
  343. /* > Invented by Sir William Rowan Hamilton */
  344. /* > Frequently used in mathematics and physics textbooks */
  345. /* - "Engineering" */
  346. /* > Widely used in aerospace engineering applications */
  347. /* SPICELIB subroutine interfaces ALWAYS use SPICE quaternions. */
  348. /* Quaternions of any other style must be converted to SPICE */
  349. /* quaternions before they are passed to SPICELIB routines. */
  350. /* Relationship between SPICE and Engineering Quaternions */
  351. /* ------------------------------------------------------ */
  352. /* Let M be a rotation matrix such that for any vector V, */
  353. /* M*V */
  354. /* is the result of rotating V by theta radians in the */
  355. /* counterclockwise direction about unit rotation axis vector A. */
  356. /* Then the SPICE quaternions representing M are */
  357. /* (+/-) ( cos(theta/2), */
  358. /* sin(theta/2) A(1), */
  359. /* sin(theta/2) A(2), */
  360. /* sin(theta/2) A(3) ) */
  361. /* while the engineering quaternions representing M are */
  362. /* (+/-) ( -sin(theta/2) A(1), */
  363. /* -sin(theta/2) A(2), */
  364. /* -sin(theta/2) A(3), */
  365. /* cos(theta/2) ) */
  366. /* For both styles of quaternions, if a quaternion q represents */
  367. /* a rotation matrix M, then -q represents M as well. */
  368. /* Given an engineering quaternion */
  369. /* QENG = ( q0, q1, q2, q3 ) */
  370. /* the equivalent SPICE quaternion is */
  371. /* QSPICE = ( q3, -q0, -q1, -q2 ) */
  372. /* Associating SPICE Quaternions with Rotation Matrices */
  373. /* ---------------------------------------------------- */
  374. /* Let FROM and TO be two right-handed reference frames, for */
  375. /* example, an inertial frame and a spacecraft-fixed frame. Let the */
  376. /* symbols */
  377. /* V , V */
  378. /* FROM TO */
  379. /* denote, respectively, an arbitrary vector expressed relative to */
  380. /* the FROM and TO frames. Let M denote the transformation matrix */
  381. /* that transforms vectors from frame FROM to frame TO; then */
  382. /* V = M * V */
  383. /* TO FROM */
  384. /* where the expression on the right hand side represents left */
  385. /* multiplication of the vector by the matrix. */
  386. /* Then if the unit-length SPICE quaternion q represents M, where */
  387. /* q = (q0, q1, q2, q3) */
  388. /* the elements of M are derived from the elements of q as follows: */
  389. /* +- -+ */
  390. /* | 2 2 | */
  391. /* | 1 - 2*( q2 + q3 ) 2*(q1*q2 - q0*q3) 2*(q1*q3 + q0*q2) | */
  392. /* | | */
  393. /* | | */
  394. /* | 2 2 | */
  395. /* M = | 2*(q1*q2 + q0*q3) 1 - 2*( q1 + q3 ) 2*(q2*q3 - q0*q1) | */
  396. /* | | */
  397. /* | | */
  398. /* | 2 2 | */
  399. /* | 2*(q1*q3 - q0*q2) 2*(q2*q3 + q0*q1) 1 - 2*( q1 + q2 ) | */
  400. /* | | */
  401. /* +- -+ */
  402. /* Note that substituting the elements of -q for those of q in the */
  403. /* right hand side leaves each element of M unchanged; this shows */
  404. /* that if a quaternion q represents a matrix M, then so does the */
  405. /* quaternion -q. */
  406. /* To map the rotation matrix M to a unit quaternion, we start by */
  407. /* decomposing the rotation matrix as a sum of symmetric */
  408. /* and skew-symmetric parts: */
  409. /* 2 */
  410. /* M = [ I + (1-cos(theta)) OMEGA ] + [ sin(theta) OMEGA ] */
  411. /* symmetric skew-symmetric */
  412. /* OMEGA is a skew-symmetric matrix of the form */
  413. /* +- -+ */
  414. /* | 0 -n3 n2 | */
  415. /* | | */
  416. /* OMEGA = | n3 0 -n1 | */
  417. /* | | */
  418. /* | -n2 n1 0 | */
  419. /* +- -+ */
  420. /* The vector N of matrix entries (n1, n2, n3) is the rotation axis */
  421. /* of M and theta is M's rotation angle. Note that N and theta */
  422. /* are not unique. */
  423. /* Let */
  424. /* C = cos(theta/2) */
  425. /* S = sin(theta/2) */
  426. /* Then the unit quaternions Q corresponding to M are */
  427. /* Q = +/- ( C, S*n1, S*n2, S*n3 ) */
  428. /* The mappings between quaternions and the corresponding rotations */
  429. /* are carried out by the SPICELIB routines */
  430. /* Q2M {quaternion to matrix} */
  431. /* M2Q {matrix to quaternion} */
  432. /* M2Q always returns a quaternion with scalar part greater than */
  433. /* or equal to zero. */
  434. /* SPICE Quaternion Multiplication Formula */
  435. /* --------------------------------------- */
  436. /* Given a SPICE quaternion */
  437. /* Q = ( q0, q1, q2, q3 ) */
  438. /* corresponding to rotation axis A and angle theta as above, we can */
  439. /* represent Q using "scalar + vector" notation as follows: */
  440. /* s = q0 = cos(theta/2) */
  441. /* v = ( q1, q2, q3 ) = sin(theta/2) * A */
  442. /* Q = s + v */
  443. /* Let Q1 and Q2 be SPICE quaternions with respective scalar */
  444. /* and vector parts s1, s2 and v1, v2: */
  445. /* Q1 = s1 + v1 */
  446. /* Q2 = s2 + v2 */
  447. /* We represent the dot product of v1 and v2 by */
  448. /* <v1, v2> */
  449. /* and the cross product of v1 and v2 by */
  450. /* v1 x v2 */
  451. /* Then the SPICE quaternion product is */
  452. /* Q1*Q2 = s1*s2 - <v1,v2> + s1*v2 + s2*v1 + (v1 x v2) */
  453. /* If Q1 and Q2 represent the rotation matrices M1 and M2 */
  454. /* respectively, then the quaternion product */
  455. /* Q1*Q2 */
  456. /* represents the matrix product */
  457. /* M1*M2 */
  458. /* $ Examples */
  459. /* Suppose that you have data packets and are prepared to produce */
  460. /* a segment of type 5 in a CK file. */
  461. /* The following code fragment could be used to add the new segment */
  462. /* to a previously opened CK file attached to HANDLE. The file must */
  463. /* have been opened with write access. */
  464. /* C */
  465. /* C Create a segment identifier. */
  466. /* C */
  467. /* SEGID = 'MY_SAMPLE_CK_TYPE_5_SEGMENT' */
  468. /* C */
  469. /* C Write the segment. */
  470. /* C */
  471. /* CALL CKW05 ( HANDLE, SUBTYP, DEGREE, BEGTIM, ENDTIM, */
  472. /* . INST, REF, AVFLAG, SEGID, N, */
  473. /* . SCLKDP, PACKTS, RATE, NINTS, STARTS ) */
  474. /* $ Restrictions */
  475. /* None. */
  476. /* $ Literature_References */
  477. /* None. */
  478. /* $ Author_and_Institution */
  479. /* N.J. Bachman (JPL) */
  480. /* W.L. Taber (JPL) */
  481. /* K.R. Gehringer (JPL) */
  482. /* J.M. Lynch (JPL) */
  483. /* $ Version */
  484. /* - SPICELIB Version 2.0.0, 08-FEB-2010 (NJB) */
  485. /* The check for non-unit quaternions has been replaced */
  486. /* with a check for zero-length quaternions. */
  487. /* - SPICELIB Version 1.1.0, 26-FEB-2008 (NJB) */
  488. /* Updated header; added information about SPICE */
  489. /* quaternion conventions. */
  490. /* Minor typo in a long error message was corrected. */
  491. /* - SPICELIB Version 1.0.1, 07-JAN-2005 (NJB) */
  492. /* Description in Detailed_Input header section of */
  493. /* constraints on BEGTIM and ENDTIM was corrected. */
  494. /* - SPICELIB Version 1.0.0, 30-AUG-2002 (NJB) (KRG) (JML) (WLT) */
  495. /* -& */
  496. /* $ Index_Entries */
  497. /* write ck type_5 data segment */
  498. /* -& */
  499. /* $ Revisions */
  500. /* - SPICELIB Version 2.0.0, 08-FEB-2010 (NJB) */
  501. /* The check for non-unit quaternions has been replaced */
  502. /* with a check for zero-length quaternions. */
  503. /* This change was made to accommodate CK generation, */
  504. /* via the non-SPICE utility MEX2KER, for European missions. */
  505. /* -& */
  506. /* SPICELIB functions */
  507. /* Local parameters */
  508. /* Packet structure parameters */
  509. /* Local variables */
  510. /* Standard SPICE error handling. */
  511. if (return_()) {
  512. return 0;
  513. } else {
  514. chkin_("CKW05", (ftnlen)5);
  515. }
  516. /* Make sure that the number of packets is positive. */
  517. if (*n < 1) {
  518. setmsg_("At least 1 packet is required for CK type 5. Number of pack"
  519. "ets supplied: #", (ftnlen)75);
  520. errint_("#", n, (ftnlen)1);
  521. sigerr_("SPICE(TOOFEWPACKETS)", (ftnlen)20);
  522. chkout_("CKW05", (ftnlen)5);
  523. return 0;
  524. }
  525. /* Make sure that there is a positive number of interpolation */
  526. /* intervals. */
  527. if (*nints <= 0) {
  528. setmsg_("# is an invalid number of interpolation intervals for type "
  529. "5.", (ftnlen)61);
  530. errint_("#", nints, (ftnlen)1);
  531. sigerr_("SPICE(INVALIDNUMINTS)", (ftnlen)21);
  532. chkout_("CKW05", (ftnlen)5);
  533. return 0;
  534. }
  535. /* Get the NAIF integer code for the reference frame. */
  536. namfrm_(ref, &refcod, ref_len);
  537. if (refcod == 0) {
  538. setmsg_("The reference frame # is not supported.", (ftnlen)39);
  539. errch_("#", ref, (ftnlen)1, ref_len);
  540. sigerr_("SPICE(INVALIDREFFRAME)", (ftnlen)22);
  541. chkout_("CKW05", (ftnlen)5);
  542. return 0;
  543. }
  544. /* Check to see if the segment identifier is too long. */
  545. if (lastnb_(segid, segid_len) > 40) {
  546. setmsg_("Segment identifier contains more than 40 characters.", (
  547. ftnlen)52);
  548. sigerr_("SPICE(SEGIDTOOLONG)", (ftnlen)19);
  549. chkout_("CKW05", (ftnlen)5);
  550. return 0;
  551. }
  552. /* Now check that all the characters in the segment identifier */
  553. /* can be printed. */
  554. i__1 = lastnb_(segid, segid_len);
  555. for (i__ = 1; i__ <= i__1; ++i__) {
  556. chrcod = *(unsigned char *)&segid[i__ - 1];
  557. if (chrcod < 32 || chrcod > 126) {
  558. setmsg_("The segment identifier contains nonprintable characters",
  559. (ftnlen)55);
  560. sigerr_("SPICE(NONPRINTABLECHARS)", (ftnlen)24);
  561. chkout_("CKW05", (ftnlen)5);
  562. return 0;
  563. }
  564. }
  565. /* Now check that the encoded SCLK times are positive and strictly */
  566. /* increasing. */
  567. /* Check that the first time is nonnegative. */
  568. if (sclkdp[0] < 0.) {
  569. setmsg_("The first SCLKDP time: # is negative.", (ftnlen)37);
  570. errdp_("#", sclkdp, (ftnlen)1);
  571. sigerr_("SPICE(INVALIDSCLKTIME)", (ftnlen)22);
  572. chkout_("CKW05", (ftnlen)5);
  573. return 0;
  574. }
  575. /* Now check that the times are ordered properly. */
  576. i__1 = *n;
  577. for (i__ = 2; i__ <= i__1; ++i__) {
  578. if (sclkdp[i__ - 1] <= sclkdp[i__ - 2]) {
  579. setmsg_("The SCLKDP times are not strictly increasing. SCLKDP(#)"
  580. " = # and SCLKDP(#) = #.", (ftnlen)78);
  581. errint_("#", &i__, (ftnlen)1);
  582. errdp_("#", &sclkdp[i__ - 1], (ftnlen)1);
  583. i__2 = i__ - 1;
  584. errint_("#", &i__2, (ftnlen)1);
  585. errdp_("#", &sclkdp[i__ - 2], (ftnlen)1);
  586. sigerr_("SPICE(TIMESOUTOFORDER)", (ftnlen)22);
  587. chkout_("CKW05", (ftnlen)5);
  588. return 0;
  589. }
  590. }
  591. /* Now check that the interval start times are ordered properly. */
  592. i__1 = *nints;
  593. for (i__ = 2; i__ <= i__1; ++i__) {
  594. if (starts[i__ - 1] <= starts[i__ - 2]) {
  595. setmsg_("The interval start times are not strictly increasing. S"
  596. "TARTS(#) = # and STARTS(#) = #.", (ftnlen)86);
  597. errint_("#", &i__, (ftnlen)1);
  598. errdp_("#", &starts[i__ - 1], (ftnlen)1);
  599. i__2 = i__ - 1;
  600. errint_("#", &i__2, (ftnlen)1);
  601. errdp_("#", &starts[i__ - 2], (ftnlen)1);
  602. sigerr_("SPICE(TIMESOUTOFORDER)", (ftnlen)22);
  603. chkout_("CKW05", (ftnlen)5);
  604. return 0;
  605. }
  606. }
  607. /* Now make sure that all of the interval start times coincide with */
  608. /* one of the times associated with the actual pointing. */
  609. i__1 = *nints;
  610. for (i__ = 1; i__ <= i__1; ++i__) {
  611. /* We know the SCLKDP array is ordered, so a binary search is */
  612. /* ok. */
  613. if (bsrchd_(&starts[i__ - 1], n, sclkdp) == 0) {
  614. setmsg_("Interval start time number # is invalid. STARTS(#) = *",
  615. (ftnlen)54);
  616. errint_("#", &i__, (ftnlen)1);
  617. errint_("#", &i__, (ftnlen)1);
  618. errdp_("*", &starts[i__ - 1], (ftnlen)1);
  619. sigerr_("SPICE(INVALIDSTARTTIME)", (ftnlen)23);
  620. chkout_("CKW05", (ftnlen)5);
  621. return 0;
  622. }
  623. }
  624. /* Set the window, packet size and angular velocity flag, all of */
  625. /* which are functions of the subtype. */
  626. if (*subtyp == 0) {
  627. winsiz = (*degree + 1) / 2;
  628. packsz = 8;
  629. } else if (*subtyp == 1) {
  630. winsiz = *degree + 1;
  631. packsz = 4;
  632. } else if (*subtyp == 2) {
  633. winsiz = (*degree + 1) / 2;
  634. packsz = 14;
  635. } else if (*subtyp == 3) {
  636. winsiz = *degree + 1;
  637. packsz = 7;
  638. } else {
  639. setmsg_("CK type 5 subtype <#> is not supported.", (ftnlen)39);
  640. errint_("#", subtyp, (ftnlen)1);
  641. sigerr_("SPICE(NOTSUPPORTED)", (ftnlen)19);
  642. chkout_("CKW05", (ftnlen)5);
  643. return 0;
  644. }
  645. /* Make sure that the quaternions are non-zero. This is just */
  646. /* a check for uninitialized data. */
  647. i__1 = *n;
  648. for (i__ = 1; i__ <= i__1; ++i__) {
  649. /* We have to address the quaternion explicitly, since the shape */
  650. /* of the packet array is not known at compile time. */
  651. addr__ = packsz * (i__ - 1) + 1;
  652. if (vzerog_(&packts[addr__ - 1], &c__4)) {
  653. setmsg_("The quaternion at index # has magnitude zero.", (ftnlen)
  654. 45);
  655. errint_("#", &i__, (ftnlen)1);
  656. sigerr_("SPICE(ZEROQUATERNION)", (ftnlen)21);
  657. chkout_("CKW05", (ftnlen)5);
  658. return 0;
  659. }
  660. }
  661. /* Make sure that the degree of the interpolating polynomials is */
  662. /* in range. */
  663. if (*degree < 1 || *degree > 15) {
  664. setmsg_("The interpolating polynomials have degree #; the valid degr"
  665. "ee range is [1, #]", (ftnlen)77);
  666. errint_("#", degree, (ftnlen)1);
  667. errint_("#", &c__15, (ftnlen)1);
  668. sigerr_("SPICE(INVALIDDEGREE)", (ftnlen)20);
  669. chkout_("CKW05", (ftnlen)5);
  670. return 0;
  671. }
  672. /* Make sure that the window size is even. If not, the input */
  673. /* DEGREE is incompatible with the subtype. */
  674. if (odd_(&winsiz)) {
  675. setmsg_("The interpolating polynomials have degree #; for CK type 5,"
  676. " the degree must be equivalent to 3 mod 4 for Hermite interp"
  677. "olation and odd for for Lagrange interpolation.", (ftnlen)166)
  678. ;
  679. errint_("#", degree, (ftnlen)1);
  680. sigerr_("SPICE(INVALIDDEGREE)", (ftnlen)20);
  681. chkout_("CKW05", (ftnlen)5);
  682. return 0;
  683. }
  684. /* If we made it this far, we're ready to start writing the segment. */
  685. /* Create the segment descriptor. */
  686. /* Assign values to the integer components of the segment descriptor. */
  687. ic[0] = *inst;
  688. ic[1] = refcod;
  689. ic[2] = 5;
  690. if (*avflag) {
  691. ic[3] = 1;
  692. } else {
  693. ic[3] = 0;
  694. }
  695. dc[0] = *begtim;
  696. dc[1] = *endtim;
  697. /* Make sure the descriptor times are in increasing order. */
  698. if (*endtim < *begtim) {
  699. setmsg_("Descriptor bounds are non-increasing: #:#", (ftnlen)41);
  700. errdp_("#", begtim, (ftnlen)1);
  701. errdp_("#", endtim, (ftnlen)1);
  702. sigerr_("SPICE(BADDESCRTIMES)", (ftnlen)20);
  703. chkout_("CKW05", (ftnlen)5);
  704. return 0;
  705. }
  706. /* Make sure that at least one time tag lies between BEGTIM and */
  707. /* ENDTIM. The first time tag not less than BEGTIM must exist */
  708. /* and must be less than or equal to ENDTIM. */
  709. i__ = lstltd_(begtim, n, sclkdp);
  710. if (i__ == *n) {
  711. setmsg_("All time tags are less than segment start time #.", (ftnlen)
  712. 49);
  713. errdp_("#", begtim, (ftnlen)1);
  714. sigerr_("SPICE(EMPTYSEGMENT)", (ftnlen)19);
  715. chkout_("CKW05", (ftnlen)5);
  716. return 0;
  717. } else if (sclkdp[i__] > *endtim) {
  718. setmsg_("No time tags lie between the segment start time # and segme"
  719. "nt end time #", (ftnlen)72);
  720. errdp_("#", begtim, (ftnlen)1);
  721. errdp_("#", endtim, (ftnlen)1);
  722. sigerr_("SPICE(EMPTYSEGMENT)", (ftnlen)19);
  723. chkout_("CKW05", (ftnlen)5);
  724. return 0;
  725. }
  726. /* The clock rate must be non-zero. */
  727. if (*rate == 0.) {
  728. setmsg_("The SCLK rate RATE was zero.", (ftnlen)28);
  729. sigerr_("SPICE(INVALIDVALUE)", (ftnlen)19);
  730. chkout_("CKW05", (ftnlen)5);
  731. return 0;
  732. }
  733. /* Now pack the segment descriptor. */
  734. dafps_(&c__2, &c__6, dc, ic, descr);
  735. /* Begin a new segment. */
  736. dafbna_(handle, descr, segid, segid_len);
  737. if (failed_()) {
  738. chkout_("CKW05", (ftnlen)5);
  739. return 0;
  740. }
  741. /* The type 5 segment structure is eloquently described by this */
  742. /* diagram from the CK Required Reading: */
  743. /* +-----------------------+ */
  744. /* | Packet 1 | */
  745. /* +-----------------------+ */
  746. /* | Packet 2 | */
  747. /* +-----------------------+ */
  748. /* . */
  749. /* . */
  750. /* . */
  751. /* +-----------------------+ */
  752. /* | Packet N | */
  753. /* +-----------------------+ */
  754. /* | Epoch 1 | */
  755. /* +-----------------------+ */
  756. /* | Epoch 2 | */
  757. /* +-----------------------+ */
  758. /* . */
  759. /* . */
  760. /* . */
  761. /* +----------------------------+ */
  762. /* | Epoch N | */
  763. /* +----------------------------+ */
  764. /* | Epoch 100 | (First directory) */
  765. /* +----------------------------+ */
  766. /* . */
  767. /* . */
  768. /* . */
  769. /* +----------------------------+ */
  770. /* | Epoch ((N-1)/100)*100 | (Last directory) */
  771. /* +----------------------------+ */
  772. /* | Start time 1 | */
  773. /* +----------------------------+ */
  774. /* | Start time 2 | */
  775. /* +----------------------------+ */
  776. /* . */
  777. /* . */
  778. /* . */
  779. /* +----------------------------+ */
  780. /* | Start time M | */
  781. /* +----------------------------+ */
  782. /* | Start time 100 | (First interval start */
  783. /* +----------------------------+ time directory) */
  784. /* . */
  785. /* . */
  786. /* . */
  787. /* +----------------------------+ */
  788. /* | Start time ((M-1)/100)*100 | (Last interval start */
  789. /* +----------------------------+ time directory) */
  790. /* | Seconds per tick | */
  791. /* +----------------------------+ */
  792. /* | Subtype code | */
  793. /* +----------------------------+ */
  794. /* | Window size | */
  795. /* +----------------------------+ */
  796. /* | Number of interp intervals | */
  797. /* +----------------------------+ */
  798. /* | Number of packets | */
  799. /* +----------------------------+ */
  800. i__1 = *n * packsz;
  801. dafada_(packts, &i__1);
  802. dafada_(sclkdp, n);
  803. i__1 = (*n - 1) / 100;
  804. for (i__ = 1; i__ <= i__1; ++i__) {
  805. dafada_(&sclkdp[i__ * 100 - 1], &c__1);
  806. }
  807. /* Now add the interval start times. */
  808. dafada_(starts, nints);
  809. /* And the directory of interval start times. The directory of */
  810. /* start times will simply be every (DIRSIZ)th start time. */
  811. i__1 = (*nints - 1) / 100;
  812. for (i__ = 1; i__ <= i__1; ++i__) {
  813. dafada_(&starts[i__ * 100 - 1], &c__1);
  814. }
  815. /* Add the SCLK rate, segment subtype, window size, interval */
  816. /* count, and packet count. */
  817. dafada_(rate, &c__1);
  818. d__1 = (doublereal) (*subtyp);
  819. dafada_(&d__1, &c__1);
  820. d__1 = (doublereal) winsiz;
  821. dafada_(&d__1, &c__1);
  822. d__1 = (doublereal) (*nints);
  823. dafada_(&d__1, &c__1);
  824. d__1 = (doublereal) (*n);
  825. dafada_(&d__1, &c__1);
  826. /* As long as nothing went wrong, end the segment. */
  827. if (! failed_()) {
  828. dafena_();
  829. }
  830. chkout_("CKW05", (ftnlen)5);
  831. return 0;
  832. } /* ckw05_ */