Chinaunix首页 | 论坛 | 博客

牛蹄印章_VX_iibull

首页　| 　博文目录　| 　关于我

iibull

  • 博客访问： 9140304
  • 博文数量： 1725
  • 博客积分： 12961
  • 博客等级： 上将
  • 技术积分： 19840
  • 用 户 组： 普通用户
  • 注册时间： 2009-01-09 11:25
个人简介

偷得浮生半桶水(半日闲), 好记性不如抄下来(烂笔头). 信息爆炸的时代, 学习是一项持续的工作.

文章分类

全部博文（1725）

文章存档

2024年（1）

2023年（26）

2022年（112）

2021年（217）

2020年（157）

2019年（192）

2018年（81）

2017年（78）

2016年（70）

2015年（52）

2014年（40）

2013年（51）

2012年（85）

2011年（45）

2010年（231）

2009年（287）

我的朋友
最近访客
推荐博文
相关博文
用加速度和角速度计算出姿态（四元数简单）

分类： Android平台

2017-10-09 19:43:04

https://wenku.baidu.com/view/8401d0e359eef8c75ebfb367.html
https://github.com/Hydr0Dr4gon/MPU9250-on-Teensy-LC/tree/master/Teensy_Tracker_Sketch


点击(此处)折叠或打开

  1. /* Dragon Labs VR Headset Tracker
  2.   by: Mark Hodgkinson
  3.   date: 20 September 2017
  4.   license: MIT Licence
  5.   Derived heavily from Kris Winer's example code for the MPU9250 and infinitellamas's code for an
  6.   Arduino Tracker for OSVR. Many thanks to both of them for making their code public and doing
  7.   most of the heavy lifting.
  8.   This is the Teensy code for the Dragon Labs VR Headset
  9.   Hardware:
  10.   Teensy LC
  11.   InvenSense MPU9250 on Sparkfun breakout
  12. */

  14. #include <i2c_t3.h>


  17. // See also MPU-9250 Register Map and Descriptions, Revision 4.0, RM-MPU-9250A-00, Rev. 1.4, 9/9/2013 for registers not listed in
  18. // above document; the MPU9250 and MPU9150 are virtually identical but the latter has a different register map
  19. //
  20. //Magnetometer Registers
  21. #define AK8963_ADDRESS 0x0C
  22. #define AK8963_WHO_AM_I 0x00 // should return 0x48
  23. #define AK8963_INFO 0x01
  24. #define AK8963_ST1 0x02 // data ready status bit 0
  25. #define AK8963_XOUT_L 0x03 // data
  26. #define AK8963_XOUT_H 0x04
  27. #define AK8963_YOUT_L 0x05
  28. #define AK8963_YOUT_H 0x06
  29. #define AK8963_ZOUT_L 0x07
  30. #define AK8963_ZOUT_H 0x08
  31. #define AK8963_ST2 0x09 // Data overflow bit 3 and data read error status bit 2
  32. #define AK8963_CNTL 0x0A // Power down (0000), single-measurement (0001), self-test (1000) and Fuse ROM (1111) modes on bits 3:0
  33. #define AK8963_ASTC 0x0C // Self test control
  34. #define AK8963_I2CDIS 0x0F // I2C disable
  35. #define AK8963_ASAX 0x10 // Fuse ROM x-axis sensitivity adjustment value
  36. #define AK8963_ASAY 0x11 // Fuse ROM y-axis sensitivity adjustment value
  37. #define AK8963_ASAZ 0x12 // Fuse ROM z-axis sensitivity adjustment value

  39. #define SELF_TEST_X_GYRO 0x00
  40. #define SELF_TEST_Y_GYRO 0x01
  41. #define SELF_TEST_Z_GYRO 0x02

  43. #define SELF_TEST_X_ACCEL 0x0D
  44. #define SELF_TEST_Y_ACCEL 0x0E
  45. #define SELF_TEST_Z_ACCEL 0x0F

  47. #define SELF_TEST_A 0x10

  49. #define XG_OFFSET_H 0x13 // User-defined trim values for gyroscope
  50. #define XG_OFFSET_L 0x14
  51. #define YG_OFFSET_H 0x15
  52. #define YG_OFFSET_L 0x16
  53. #define ZG_OFFSET_H 0x17
  54. #define ZG_OFFSET_L 0x18
  55. #define SMPLRT_DIV 0x19
  56. #define CONFIG 0x1A
  57. #define GYRO_CONFIG 0x1B
  58. #define ACCEL_CONFIG 0x1C
  59. #define ACCEL_CONFIG2 0x1D
  60. #define LP_ACCEL_ODR 0x1E
  61. #define WOM_THR 0x1F

  63. #define MOT_DUR 0x20 // Duration counter threshold for motion interrupt generation, 1 kHz rate, LSB = 1 ms
  64. #define ZMOT_THR 0x21 // Zero-motion detection threshold bits [7:0]
  65. #define ZRMOT_DUR 0x22 // Duration counter threshold for zero motion interrupt generation, 16 Hz rate, LSB = 64 ms

  67. #define FIFO_EN 0x23
  68. #define I2C_MST_CTRL 0x24
  69. #define I2C_SLV0_ADDR 0x25
  70. #define I2C_SLV0_REG 0x26
  71. #define I2C_SLV0_CTRL 0x27
  72. #define I2C_SLV1_ADDR 0x28
  73. #define I2C_SLV1_REG 0x29
  74. #define I2C_SLV1_CTRL 0x2A
  75. #define I2C_SLV2_ADDR 0x2B
  76. #define I2C_SLV2_REG 0x2C
  77. #define I2C_SLV2_CTRL 0x2D
  78. #define I2C_SLV3_ADDR 0x2E
  79. #define I2C_SLV3_REG 0x2F
  80. #define I2C_SLV3_CTRL 0x30
  81. #define I2C_SLV4_ADDR 0x31
  82. #define I2C_SLV4_REG 0x32
  83. #define I2C_SLV4_DO 0x33
  84. #define I2C_SLV4_CTRL 0x34
  85. #define I2C_SLV4_DI 0x35
  86. #define I2C_MST_STATUS 0x36
  87. #define INT_PIN_CFG 0x37
  88. #define INT_ENABLE 0x38
  89. #define DMP_INT_STATUS 0x39 // Check DMP interrupt
  90. #define INT_STATUS 0x3A
  91. #define ACCEL_XOUT_H 0x3B
  92. #define ACCEL_XOUT_L 0x3C
  93. #define ACCEL_YOUT_H 0x3D
  94. #define ACCEL_YOUT_L 0x3E
  95. #define ACCEL_ZOUT_H 0x3F
  96. #define ACCEL_ZOUT_L 0x40
  97. #define TEMP_OUT_H 0x41
  98. #define TEMP_OUT_L 0x42
  99. #define GYRO_XOUT_H 0x43
  100. #define GYRO_XOUT_L 0x44
  101. #define GYRO_YOUT_H 0x45
  102. #define GYRO_YOUT_L 0x46
  103. #define GYRO_ZOUT_H 0x47
  104. #define GYRO_ZOUT_L 0x48
  105. #define EXT_SENS_DATA_00 0x49
  106. #define EXT_SENS_DATA_01 0x4A
  107. #define EXT_SENS_DATA_02 0x4B
  108. #define EXT_SENS_DATA_03 0x4C
  109. #define EXT_SENS_DATA_04 0x4D
  110. #define EXT_SENS_DATA_05 0x4E
  111. #define EXT_SENS_DATA_06 0x4F
  112. #define EXT_SENS_DATA_07 0x50
  113. #define EXT_SENS_DATA_08 0x51
  114. #define EXT_SENS_DATA_09 0x52
  115. #define EXT_SENS_DATA_10 0x53
  116. #define EXT_SENS_DATA_11 0x54
  117. #define EXT_SENS_DATA_12 0x55
  118. #define EXT_SENS_DATA_13 0x56
  119. #define EXT_SENS_DATA_14 0x57
  120. #define EXT_SENS_DATA_15 0x58
  121. #define EXT_SENS_DATA_16 0x59
  122. #define EXT_SENS_DATA_17 0x5A
  123. #define EXT_SENS_DATA_18 0x5B
  124. #define EXT_SENS_DATA_19 0x5C
  125. #define EXT_SENS_DATA_20 0x5D
  126. #define EXT_SENS_DATA_21 0x5E
  127. #define EXT_SENS_DATA_22 0x5F
  128. #define EXT_SENS_DATA_23 0x60
  129. #define MOT_DETECT_STATUS 0x61
  130. #define I2C_SLV0_DO 0x63
  131. #define I2C_SLV1_DO 0x64
  132. #define I2C_SLV2_DO 0x65
  133. #define I2C_SLV3_DO 0x66
  134. #define I2C_MST_DELAY_CTRL 0x67
  135. #define SIGNAL_PATH_RESET 0x68
  136. #define MOT_DETECT_CTRL 0x69
  137. #define USER_CTRL 0x6A // Bit 7 enable DMP, bit 3 reset DMP
  138. #define PWR_MGMT_1 0x6B // Device defaults to the SLEEP mode
  139. #define PWR_MGMT_2 0x6C
  140. #define DMP_BANK 0x6D // Activates a specific bank in the DMP
  141. #define DMP_RW_PNT 0x6E // Set read/write pointer to a specific start address in specified DMP bank
  142. #define DMP_REG 0x6F // Register in DMP from which to read or to which to write
  143. #define DMP_REG_1 0x70
  144. #define DMP_REG_2 0x71
  145. #define FIFO_COUNTH 0x72
  146. #define FIFO_COUNTL 0x73
  147. #define FIFO_R_W 0x74
  148. #define WHO_AM_I_MPU9250 0x75 // Should return 0x71
  149. #define XA_OFFSET_H 0x77
  150. #define XA_OFFSET_L 0x78
  151. #define YA_OFFSET_H 0x7A
  152. #define YA_OFFSET_L 0x7B
  153. #define ZA_OFFSET_H 0x7D
  154. #define ZA_OFFSET_L 0x7E

  156. // Using the MSENSR-9250 breakout board, ADO is set to 0
  157. // Seven-bit device address is 110100 for ADO = 0 and 110101 for ADO = 1
  158. //#define ADO 0
  159. #if ADO
  160. #define MPU9250_ADDRESS 0x69 // Device address when ADO = 1
  161. #else
  162. #define MPU9250_ADDRESS 0x68 // Device address when ADO = 0
  163. #define AK8963_ADDRESS 0x0C // Address of magnetometer
  164. #endif

  166. // Set initial input parameters
  167. enum Ascale {
  168.   AFS_2G = 0,
  169.   AFS_4G,
  170.   AFS_8G,
  171.   AFS_16G
  172. };

  174. enum Gscale {
  175.   GFS_250DPS = 0,
  176.   GFS_500DPS,
  177.   GFS_1000DPS,
  178.   GFS_2000DPS
  179. };

  181. enum Mscale {
  182.   MFS_14BITS = 0, // 0.6 mG per LSB
  183.   MFS_16BITS // 0.15 mG per LSB
  184. };

  186. // Specify sensor full scale
  187. uint8_t Gscale = GFS_250DPS;
  188. uint8_t Ascale = AFS_2G;
  189. uint8_t Mscale = MFS_16BITS; // Choose either 14-bit or 16-bit magnetometer resolution
  190. uint8_t Mmode = 0x02; // 2 for 8 Hz, 6 for 100 Hz continuous magnetometer data read
  191. float aRes, gRes, mRes; // scale resolutions per LSB for the sensors

  193. //// Pin definitions
  194. //int intPin = 12; // These can be changed, 2 and 3 are the Arduinos ext int pins
  195. //int adoPin = 8;
  196. //int LEDPin = 13;

  198. int16_t accelCount[3]; // Stores the 16-bit signed accelerometer sensor output
  199. int16_t gyroCount[3]; // Stores the 16-bit signed gyro sensor output
  200. int16_t magCount[3]; // Stores the 16-bit signed magnetometer sensor output
  201. float magCalibration[3] = {0, 0, 0}, magbias[3] = {0, 0, 0}; // Factory mag calibration and mag bias
  202. float gyroBias[3] = {0, 0, 0}, accelBias[3] = {0, 0, 0}; // Bias corrections for gyro and accelerometer
  203. int16_t tempCount; // temperature raw count output
  204. float temperature; // Stores the real internal chip temperature in degrees Celsius
  205. float SelfTest[6]; // holds results of gyro and accelerometer self test

  207. // global constants for 9 DoF fusion and AHRS (Attitude and Heading Reference System)
  208. float GyroMeasError = PI * (40.0f / 180.0f); // gyroscope measurement error in rads/s (start at 40 deg/s)
  209. float GyroMeasDrift = PI * (0.0f / 180.0f); // gyroscope measurement drift in rad/s/s (start at 0.0 deg/s/s)

  211. // There is a tradeoff in the beta parameter between accuracy and response speed.
  212. // In the original Madgwick study, beta of 0.041 (corresponding to GyroMeasError of 2.7 degrees/s) was found to give optimal accuracy.
  213. // However, with this value, the LSM9SD0 response time is about 10 seconds to a stable initial quaternion.
  214. // Subsequent changes also require a longish lag time to a stable output, not fast enough for a quadcopter or robot
  215. // By increasing beta (GyroMeasError) by about a factor of fifteen, the response time constant is reduced to ~2 sec
  216. // I haven't noticed any reduction in solution accuracy. This is essentially the I coefficient in a PID control sense;
  217. // the bigger the feedback coefficient, the faster the solution converges, usually at the expense of accuracy.
  218. // In any case, this is the free parameter in the Madgwick filtering and fusion scheme.
  219. float beta = sqrt(3.0f / 4.0f) * GyroMeasError; // compute beta
  220. float zeta = sqrt(3.0f / 4.0f) * GyroMeasDrift; // compute zeta, the other free parameter in the Madgwick scheme usually set to a small or zero value
  221. #define Kp 2.0f * 5.0f // these are the free parameters in the Mahony filter and fusion scheme, Kp for proportional feedback, Ki for integral
  222. #define Ki 0.0f

  224. uint32_t delt_t = 0; // used to control display output rate
  225. uint32_t count = 0, sumCount = 0; // used to control display output rate
  226. //float pitch, yaw, roll;
  227. float deltat = 0.0f, sum = 0.0f; // integration interval for both filter schemes
  228. uint32_t lastUpdate = 0, firstUpdate = 0; // used to calculate integration interval
  229. uint32_t Now = 0; // used to calculate integration interval

  231. float ax, ay, az, gx, gy, gz, mx, my, mz; // variables to hold latest sensor data values
  232. float q[4] = {1.0f, 0.0f, 0.0f, 0.0f}; // vector to hold quaternion
  233. float eInt[3] = {0.0f, 0.0f, 0.0f}; // vector to hold integral error for Mahony method


  236. void setup()
  237. {
  238.   // Setup for Master mode, pins 18/19, external pullups, 400kHz
  239.   Wire1.begin(I2C_MASTER, 0x68, I2C_PINS_22_23, I2C_PULLUP_EXT, I2C_RATE_100);
  240.   Serial.begin(115200);

  242.   // Read the WHO_AM_I register, this is a good test of communication
  243.   byte c = readByte(MPU9250_ADDRESS, WHO_AM_I_MPU9250); // Read WHO_AM_I register for MPU-9250
  244.   //Serial.print("MPU9250 "); Serial.print("I AM "); Serial.print(c, HEX); Serial.print(" I should be "); Serial.println(0x71, HEX);
  245.   delay(5000);

  247.   if (c == 0x71) // WHO_AM_I should always be 0x68
  248.   {
  249.     Serial.println("C:Begin"); //Signal calibrating sensors

  251.     MPU9250SelfTest(SelfTest); // Start by performing self test and reporting values
  252.     delay(5000);

  254.     calibrateMPU9250(gyroBias, accelBias); // Calibrate gyro and accelerometers, load biases in bias registers
  255.     delay(1000);

  257.     initMPU9250(); // Initialize device for active mode read of acclerometer, gyroscope, and temperature

  259.     delay(1000);

  261.     // Get magnetometer calibration from AK8963 ROM
  262.     initAK8963(magCalibration); // Initialize device for active mode read of magnetometer
  263.     
  264.     Serial.println("C:End"); //Signals done with calibration
  265.     
  266.   }
  267.   else
  268.   {
  269.     Serial.print("E:");
  270.     Serial.print("Could not connect to MPU9250: 0x");
  271.     Serial.println(c, HEX);
  272.     while (1) ; // Loop forever if communication doesn't happen
  273.   }
  274. }

  276. void loop()
  277. {
  278.   // If intPin goes high, all data registers have new data
  279.   if (readByte(MPU9250_ADDRESS, INT_STATUS) & 0x01) { // On interrupt, check if data ready interrupt
  280.     readAccelData(accelCount); // Read the x/y/z adc values
  281.     getAres();

  283.     // Now we'll calculate the accleration value into actual g's
  284.     ax = (float)accelCount[0] * aRes; // - accelBias[0]; // get actual g value, this depends on scale being set
  285.     ay = (float)accelCount[1] * aRes; // - accelBias[1];
  286.     az = (float)accelCount[2] * aRes; // - accelBias[2];

  288.     readGyroData(gyroCount); // Read the x/y/z adc values
  289.     getGres();

  291.     // Calculate the gyro value into actual degrees per second
  292.     gx = (float)gyroCount[0] * gRes; // get actual gyro value, this depends on scale being set
  293.     gy = (float)gyroCount[1] * gRes;
  294.     gz = (float)gyroCount[2] * gRes;

  296.     readMagData(magCount); // Read the x/y/z adc values
  297.     getMres();
  298.     magbias[0] = +470.; // User environmental x-axis correction in milliGauss, should be automatically calculated
  299.     magbias[1] = +120.; // User environmental x-axis correction in milliGauss
  300.     magbias[2] = +125.; // User environmental x-axis correction in milliGauss

  302.     // Calculate the magnetometer values in milliGauss
  303.     // Include factory calibration per data sheet and user environmental corrections
  304.     mx = (float)magCount[0] * mRes * magCalibration[0] - magbias[0]; // get actual magnetometer value, this depends on scale being set
  305.     my = (float)magCount[1] * mRes * magCalibration[1] - magbias[1];
  306.     mz = (float)magCount[2] * mRes * magCalibration[2] - magbias[2];
  307.   }

  309.   Now = micros();
  310.   deltat = ((Now - lastUpdate) / 1000000.0f); // set integration time by time elapsed since last filter update
  311.   lastUpdate = Now;

  313.   sum += deltat; // sum for averaging filter update rate
  314.   sumCount++;

  316.   // Sensors x (y)-axis of the accelerometer is aligned with the y (x)-axis of the magnetometer;
  317.   // the magnetometer z-axis (+ down) is opposite to z-axis (+ up) of accelerometer and
  318.   // We have to make some allowance for this orientationmismatch in feeding the output to the quaternion filter.
  319.   
  320.   MadgwickQuaternionUpdate(ax, ay, az, gx * PI / 180.0f, gy * PI / 180.0f, gz * PI / 180.0f, my, mx, mz);

  322.       Serial.print("Q:");
  323.       Serial.print(q[0], 5); //Print qw
  324.       Serial.print(",");
  325.       Serial.print(q[1], 5); //Print qx
  326.       Serial.print(",");
  327.       Serial.print(q[2], 5); //Print qy
  328.       Serial.print(",");
  329.       Serial.println(q[3], 5); //Print qz
  330. }



  334. //===================================================================================================================
  335. //====== Set of useful function to access acceleration. gyroscope, magnetometer, and temperature data
  336. //===================================================================================================================

  338. void getMres() {
  339.   switch (Mscale)
  340.   {
  341.     // Possible magnetometer scales (and their register bit settings) are:
  342.     // 14 bit resolution (0) and 16 bit resolution (1)
  343.     case MFS_14BITS:
  344.       mRes = 10.*4912. / 8190.; // Proper scale to return milliGauss
  345.       break;
  346.     case MFS_16BITS:
  347.       mRes = 10.*4912. / 32760.0; // Proper scale to return milliGauss
  348.       break;
  349.   }
  350. }

  352. void getGres() {
  353.   switch (Gscale)
  354.   {
  355.     // Possible gyro scales (and their register bit settings) are:
  356.     // 250 DPS (00), 500 DPS (01), 1000 DPS (10), and 2000 DPS (11).
  357.     // Here's a bit of an algorith to calculate DPS/(ADC tick) based on that 2-bit value:
  358.     case GFS_250DPS:
  359.       gRes = 250.0 / 32768.0;
  360.       break;
  361.     case GFS_500DPS:
  362.       gRes = 500.0 / 32768.0;
  363.       break;
  364.     case GFS_1000DPS:
  365.       gRes = 1000.0 / 32768.0;
  366.       break;
  367.     case GFS_2000DPS:
  368.       gRes = 2000.0 / 32768.0;
  369.       break;
  370.   }
  371. }

  373. void getAres() {
  374.   switch (Ascale)
  375.   {
  376.     // Possible accelerometer scales (and their register bit settings) are:
  377.     // 2 Gs (00), 4 Gs (01), 8 Gs (10), and 16 Gs (11).
  378.     // Here's a bit of an algorith to calculate DPS/(ADC tick) based on that 2-bit value:
  379.     case AFS_2G:
  380.       aRes = 2.0 / 32768.0;
  381.       break;
  382.     case AFS_4G:
  383.       aRes = 4.0 / 32768.0;
  384.       break;
  385.     case AFS_8G:
  386.       aRes = 8.0 / 32768.0;
  387.       break;
  388.     case AFS_16G:
  389.       aRes = 16.0 / 32768.0;
  390.       break;
  391.   }
  392. }


  395. void readAccelData(int16_t * destination)
  396. {
  397.   uint8_t rawData[6]; // x/y/z accel register data stored here
  398.   readBytes(MPU9250_ADDRESS, ACCEL_XOUT_H, 6, &rawData[0]); // Read the six raw data registers into data array
  399.   destination[0] = ((int16_t)rawData[0] << 8) | rawData[1] ; // Turn the MSB and LSB into a signed 16-bit value
  400.   destination[1] = ((int16_t)rawData[2] << 8) | rawData[3] ;
  401.   destination[2] = ((int16_t)rawData[4] << 8) | rawData[5] ;
  402. }


  405. void readGyroData(int16_t * destination)
  406. {
  407.   uint8_t rawData[6]; // x/y/z gyro register data stored here
  408.   readBytes(MPU9250_ADDRESS, GYRO_XOUT_H, 6, &rawData[0]); // Read the six raw data registers sequentially into data array
  409.   destination[0] = ((int16_t)rawData[0] << 8) | rawData[1] ; // Turn the MSB and LSB into a signed 16-bit value
  410.   destination[1] = ((int16_t)rawData[2] << 8) | rawData[3] ;
  411.   destination[2] = ((int16_t)rawData[4] << 8) | rawData[5] ;
  412. }

  414. void readMagData(int16_t * destination)
  415. {
  416.   uint8_t rawData[7]; // x/y/z gyro register data, ST2 register stored here, must read ST2 at end of data acquisition
  417.   if (readByte(AK8963_ADDRESS, AK8963_ST1) & 0x01) { // wait for magnetometer data ready bit to be set
  418.     readBytes(AK8963_ADDRESS, AK8963_XOUT_L, 7, &rawData[0]); // Read the six raw data and ST2 registers sequentially into data array
  419.     uint8_t c = rawData[6]; // End data read by reading ST2 register
  420.     if (!(c & 0x08)) { // Check if magnetic sensor overflow set, if not then report data
  421.       destination[0] = ((int16_t)rawData[1] << 8) | rawData[0] ; // Turn the MSB and LSB into a signed 16-bit value
  422.       destination[1] = ((int16_t)rawData[3] << 8) | rawData[2] ; // Data stored as little Endian
  423.       destination[2] = ((int16_t)rawData[5] << 8) | rawData[4] ;
  424.     }
  425.   }
  426. }

  428. void initAK8963(float * destination)
  429. {
  430.   // First extract the factory calibration for each magnetometer axis
  431.   uint8_t rawData[3]; // x/y/z gyro calibration data stored here
  432.   writeByte(AK8963_ADDRESS, AK8963_CNTL, 0x00); // Power down magnetometer
  433.   delay(10);
  434.   writeByte(AK8963_ADDRESS, AK8963_CNTL, 0x0F); // Enter Fuse ROM access mode
  435.   delay(10);
  436.   readBytes(AK8963_ADDRESS, AK8963_ASAX, 3, &rawData[0]); // Read the x-, y-, and z-axis calibration values
  437.   destination[0] = (float)(rawData[0] - 128) / 256. + 1.; // Return x-axis sensitivity adjustment values, etc.
  438.   destination[1] = (float)(rawData[1] - 128) / 256. + 1.;
  439.   destination[2] = (float)(rawData[2] - 128) / 256. + 1.;
  440.   writeByte(AK8963_ADDRESS, AK8963_CNTL, 0x00); // Power down magnetometer
  441.   delay(10);
  442.   // Configure the magnetometer for continuous read and highest resolution
  443.   // set Mscale bit 4 to 1 (0) to enable 16 (14) bit resolution in CNTL register,
  444.   // and enable continuous mode data acquisition Mmode (bits [3:0]), 0010 for 8 Hz and 0110 for 100 Hz sample rates
  445.   writeByte(AK8963_ADDRESS, AK8963_CNTL, Mscale << 4 | Mmode); // Set magnetometer data resolution and sample ODR
  446.   delay(10);
  447. }


  450. void initMPU9250()
  451. {
  452.   // wake up device
  453.   writeByte(MPU9250_ADDRESS, PWR_MGMT_1, 0x00); // Clear sleep mode bit (6), enable all sensors
  454.   delay(100); // Wait for all registers to reset

  456.   // get stable time source
  457.   writeByte(MPU9250_ADDRESS, PWR_MGMT_1, 0x01); // Auto select clock source to be PLL gyroscope reference if ready else
  458.   delay(200);

  460.   // Configure Gyro and Thermometer
  461.   // Disable FSYNC and set thermometer and gyro bandwidth to 41 and 42 Hz, respectively;
  462.   // minimum delay time for this setting is 5.9 ms, which means sensor fusion update rates cannot
  463.   // be higher than 1 / 0.0059 = 170 Hz
  464.   // DLPF_CFG = bits 2:0 = 011; this limits the sample rate to 1000 Hz for both
  465.   // With the MPU9250, it is possible to get gyro sample rates of 32 kHz (!), 8 kHz, or 1 kHz
  466.   writeByte(MPU9250_ADDRESS, CONFIG, 0x03);

  468.   // Set sample rate = gyroscope output rate/(1 + SMPLRT_DIV)
  469.   writeByte(MPU9250_ADDRESS, SMPLRT_DIV, 0x04); // Use a 200 Hz rate; a rate consistent with the filter update rate
  470.   // determined inset in CONFIG above

  472.   // Set gyroscope full scale range
  473.   // Range selects FS_SEL and GFS_SEL are 0 - 3, so 2-bit values are left-shifted into positions 4:3
  474.   uint8_t c = readByte(MPU9250_ADDRESS, GYRO_CONFIG); // get current GYRO_CONFIG register value
  475.   // c = c & ~0xE0; // Clear self-test bits [7:5]
  476.   c = c & ~0x03; // Clear Fchoice bits [1:0]
  477.   c = c & ~0x18; // Clear GFS bits [4:3]
  478.   c = c | Gscale << 3; // Set full scale range for the gyro
  479.   // c =| 0x00; // Set Fchoice for the gyro to 11 by writing its inverse to bits 1:0 of GYRO_CONFIG
  480.   writeByte(MPU9250_ADDRESS, GYRO_CONFIG, c ); // Write new GYRO_CONFIG value to register

  482.   // Set accelerometer full-scale range configuration
  483.   c = readByte(MPU9250_ADDRESS, ACCEL_CONFIG); // get current ACCEL_CONFIG register value
  484.   // c = c & ~0xE0; // Clear self-test bits [7:5]
  485.   c = c & ~0x18; // Clear AFS bits [4:3]
  486.   c = c | Ascale << 3; // Set full scale range for the accelerometer
  487.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG, c); // Write new ACCEL_CONFIG register value

  489.   // Set accelerometer sample rate configuration
  490.   // It is possible to get a 4 kHz sample rate from the accelerometer by choosing 1 for
  491.   // accel_fchoice_b bit [3]; in this case the bandwidth is 1.13 kHz
  492.   c = readByte(MPU9250_ADDRESS, ACCEL_CONFIG2); // get current ACCEL_CONFIG2 register value
  493.   c = c & ~0x0F; // Clear accel_fchoice_b (bit 3) and A_DLPFG (bits [2:0])
  494.   c = c | 0x03; // Set accelerometer rate to 1 kHz and bandwidth to 41 Hz
  495.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG2, c); // Write new ACCEL_CONFIG2 register value

  497.   // The accelerometer, gyro, and thermometer are set to 1 kHz sample rates,
  498.   // but all these rates are further reduced by a factor of 5 to 200 Hz because of the SMPLRT_DIV setting

  500.   // Configure Interrupts and Bypass Enable
  501.   // Set interrupt pin active high, push-pull, hold interrupt pin level HIGH until interrupt cleared,
  502.   // clear on read of INT_STATUS, and enable I2C_BYPASS_EN so additional chips
  503.   // can join the I2C bus and all can be controlled by the Arduino as master
  504.   writeByte(MPU9250_ADDRESS, INT_PIN_CFG, 0x22);
  505.   writeByte(MPU9250_ADDRESS, INT_ENABLE, 0x01); // Enable data ready (bit 0) interrupt
  506.   delay(100);
  507. }


  510. // Function which accumulates gyro and accelerometer data after device initialization. It calculates the average
  511. // of the at-rest readings and then loads the resulting offsets into accelerometer and gyro bias registers.
  512. void calibrateMPU9250(float * dest1, float * dest2)
  513. {
  514.   uint8_t data[12]; // data array to hold accelerometer and gyro x, y, z, data
  515.   uint16_t ii, packet_count, fifo_count;
  516.   int32_t gyro_bias[3] = {0, 0, 0}, accel_bias[3] = {0, 0, 0};

  518.   // reset device
  519.   writeByte(MPU9250_ADDRESS, PWR_MGMT_1, 0x80); // Write a one to bit 7 reset bit; toggle reset device
  520.   delay(100);

  522.   // get stable time source; Auto select clock source to be PLL gyroscope reference if ready
  523.   // else use the internal oscillator, bits 2:0 = 001
  524.   writeByte(MPU9250_ADDRESS, PWR_MGMT_1, 0x01);
  525.   writeByte(MPU9250_ADDRESS, PWR_MGMT_2, 0x00);
  526.   delay(200);

  528.   // Configure device for bias calculation
  529.   writeByte(MPU9250_ADDRESS, INT_ENABLE, 0x00); // Disable all interrupts
  530.   writeByte(MPU9250_ADDRESS, FIFO_EN, 0x00); // Disable FIFO
  531.   writeByte(MPU9250_ADDRESS, PWR_MGMT_1, 0x00); // Turn on internal clock source
  532.   writeByte(MPU9250_ADDRESS, I2C_MST_CTRL, 0x00); // Disable I2C master
  533.   writeByte(MPU9250_ADDRESS, USER_CTRL, 0x00); // Disable FIFO and I2C master modes
  534.   writeByte(MPU9250_ADDRESS, USER_CTRL, 0x0C); // Reset FIFO and DMP
  535.   delay(15);

  537.   // Configure MPU6050 gyro and accelerometer for bias calculation
  538.   writeByte(MPU9250_ADDRESS, CONFIG, 0x01); // Set low-pass filter to 188 Hz
  539.   writeByte(MPU9250_ADDRESS, SMPLRT_DIV, 0x00); // Set sample rate to 1 kHz
  540.   writeByte(MPU9250_ADDRESS, GYRO_CONFIG, 0x00); // Set gyro full-scale to 250 degrees per second, maximum sensitivity
  541.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG, 0x00); // Set accelerometer full-scale to 2 g, maximum sensitivity

  543.   uint16_t gyrosensitivity = 131; // = 131 LSB/degrees/sec
  544.   uint16_t accelsensitivity = 16384; // = 16384 LSB/g

  546.   // Configure FIFO to capture accelerometer and gyro data for bias calculation
  547.   writeByte(MPU9250_ADDRESS, USER_CTRL, 0x40); // Enable FIFO
  548.   writeByte(MPU9250_ADDRESS, FIFO_EN, 0x78); // Enable gyro and accelerometer sensors for FIFO (max size 512 bytes in MPU-9150)
  549.   delay(40); // accumulate 40 samples in 40 milliseconds = 480 bytes

  551.   // At end of sample accumulation, turn off FIFO sensor read
  552.   writeByte(MPU9250_ADDRESS, FIFO_EN, 0x00); // Disable gyro and accelerometer sensors for FIFO
  553.   readBytes(MPU9250_ADDRESS, FIFO_COUNTH, 2, &data[0]); // read FIFO sample count
  554.   fifo_count = ((uint16_t)data[0] << 8) | data[1];
  555.   packet_count = fifo_count / 12; // How many sets of full gyro and accelerometer data for averaging

  557.   for (ii = 0; ii < packet_count; ii++) {
  558.     int16_t accel_temp[3] = {0, 0, 0}, gyro_temp[3] = {0, 0, 0};
  559.     readBytes(MPU9250_ADDRESS, FIFO_R_W, 12, &data[0]); // read data for averaging
  560.     accel_temp[0] = (int16_t) (((int16_t)data[0] << 8) | data[1] ) ; // Form signed 16-bit integer for each sample in FIFO
  561.     accel_temp[1] = (int16_t) (((int16_t)data[2] << 8) | data[3] ) ;
  562.     accel_temp[2] = (int16_t) (((int16_t)data[4] << 8) | data[5] ) ;
  563.     gyro_temp[0] = (int16_t) (((int16_t)data[6] << 8) | data[7] ) ;
  564.     gyro_temp[1] = (int16_t) (((int16_t)data[8] << 8) | data[9] ) ;
  565.     gyro_temp[2] = (int16_t) (((int16_t)data[10] << 8) | data[11]) ;

  567.     accel_bias[0] += (int32_t) accel_temp[0]; // Sum individual signed 16-bit biases to get accumulated signed 32-bit biases
  568.     accel_bias[1] += (int32_t) accel_temp[1];
  569.     accel_bias[2] += (int32_t) accel_temp[2];
  570.     gyro_bias[0] += (int32_t) gyro_temp[0];
  571.     gyro_bias[1] += (int32_t) gyro_temp[1];
  572.     gyro_bias[2] += (int32_t) gyro_temp[2];

  574.   }
  575.   accel_bias[0] /= (int32_t) packet_count; // Normalize sums to get average count biases
  576.   accel_bias[1] /= (int32_t) packet_count;
  577.   accel_bias[2] /= (int32_t) packet_count;
  578.   gyro_bias[0] /= (int32_t) packet_count;
  579.   gyro_bias[1] /= (int32_t) packet_count;
  580.   gyro_bias[2] /= (int32_t) packet_count;

  582.   if (accel_bias[2] > 0L) {
  583.     accel_bias[2] -= (int32_t) accelsensitivity; // Remove gravity from the z-axis accelerometer bias calculation
  584.   }
  585.   else {
  586.     accel_bias[2] += (int32_t) accelsensitivity;
  587.   }

  589.   // Construct the gyro biases for push to the hardware gyro bias registers, which are reset to zero upon device startup
  590.   data[0] = (-gyro_bias[0] / 4 >> 8) & 0xFF; // Divide by 4 to get 32.9 LSB per deg/s to conform to expected bias input format
  591.   data[1] = (-gyro_bias[0] / 4) & 0xFF; // Biases are additive, so change sign on calculated average gyro biases
  592.   data[2] = (-gyro_bias[1] / 4 >> 8) & 0xFF;
  593.   data[3] = (-gyro_bias[1] / 4) & 0xFF;
  594.   data[4] = (-gyro_bias[2] / 4 >> 8) & 0xFF;
  595.   data[5] = (-gyro_bias[2] / 4) & 0xFF;

  597.   // Push gyro biases to hardware registers
  598.   writeByte(MPU9250_ADDRESS, XG_OFFSET_H, data[0]);
  599.   writeByte(MPU9250_ADDRESS, XG_OFFSET_L, data[1]);
  600.   writeByte(MPU9250_ADDRESS, YG_OFFSET_H, data[2]);
  601.   writeByte(MPU9250_ADDRESS, YG_OFFSET_L, data[3]);
  602.   writeByte(MPU9250_ADDRESS, ZG_OFFSET_H, data[4]);
  603.   writeByte(MPU9250_ADDRESS, ZG_OFFSET_L, data[5]);

  605.   // Output scaled gyro biases for display in the main program
  606.   dest1[0] = (float) gyro_bias[0] / (float) gyrosensitivity;
  607.   dest1[1] = (float) gyro_bias[1] / (float) gyrosensitivity;
  608.   dest1[2] = (float) gyro_bias[2] / (float) gyrosensitivity;

  610.   // Construct the accelerometer biases for push to the hardware accelerometer bias registers. These registers contain
  611.   // factory trim values which must be added to the calculated accelerometer biases; on boot up these registers will hold
  612.   // non-zero values. In addition, bit 0 of the lower byte must be preserved since it is used for temperature
  613.   // compensation calculations. Accelerometer bias registers expect bias input as 2048 LSB per g, so that
  614.   // the accelerometer biases calculated above must be divided by 8.

  616.   int32_t accel_bias_reg[3] = {0, 0, 0}; // A place to hold the factory accelerometer trim biases
  617.   readBytes(MPU9250_ADDRESS, XA_OFFSET_H, 2, &data[0]); // Read factory accelerometer trim values
  618.   accel_bias_reg[0] = (int32_t) (((int16_t)data[0] << 8) | data[1]);
  619.   readBytes(MPU9250_ADDRESS, YA_OFFSET_H, 2, &data[0]);
  620.   accel_bias_reg[1] = (int32_t) (((int16_t)data[0] << 8) | data[1]);
  621.   readBytes(MPU9250_ADDRESS, ZA_OFFSET_H, 2, &data[0]);
  622.   accel_bias_reg[2] = (int32_t) (((int16_t)data[0] << 8) | data[1]);

  624.   uint32_t mask = 1uL; // Define mask for temperature compensation bit 0 of lower byte of accelerometer bias registers
  625.   uint8_t mask_bit[3] = {0, 0, 0}; // Define array to hold mask bit for each accelerometer bias axis

  627.   for (ii = 0; ii < 3; ii++) {
  628.     if ((accel_bias_reg[ii] & mask)) mask_bit[ii] = 0x01; // If temperature compensation bit is set, record that fact in mask_bit
  629.   }

  631.   // Construct total accelerometer bias, including calculated average accelerometer bias from above
  632.   accel_bias_reg[0] -= (accel_bias[0] / 8); // Subtract calculated averaged accelerometer bias scaled to 2048 LSB/g (16 g full scale)
  633.   accel_bias_reg[1] -= (accel_bias[1] / 8);
  634.   accel_bias_reg[2] -= (accel_bias[2] / 8);

  636.   data[0] = (accel_bias_reg[0] >> 8) & 0xFF;
  637.   data[1] = (accel_bias_reg[0]) & 0xFF;
  638.   data[1] = data[1] | mask_bit[0]; // preserve temperature compensation bit when writing back to accelerometer bias registers
  639.   data[2] = (accel_bias_reg[1] >> 8) & 0xFF;
  640.   data[3] = (accel_bias_reg[1]) & 0xFF;
  641.   data[3] = data[3] | mask_bit[1]; // preserve temperature compensation bit when writing back to accelerometer bias registers
  642.   data[4] = (accel_bias_reg[2] >> 8) & 0xFF;
  643.   data[5] = (accel_bias_reg[2]) & 0xFF;
  644.   data[5] = data[5] | mask_bit[2]; // preserve temperature compensation bit when writing back to accelerometer bias registers

  646.   // Apparently this is not working for the acceleration biases in the MPU-9250
  647.   // Are we handling the temperature correction bit properly?
  648.   // Push accelerometer biases to hardware registers
  649.   writeByte(MPU9250_ADDRESS, XA_OFFSET_H, data[0]);
  650.   writeByte(MPU9250_ADDRESS, XA_OFFSET_L, data[1]);
  651.   writeByte(MPU9250_ADDRESS, YA_OFFSET_H, data[2]);
  652.   writeByte(MPU9250_ADDRESS, YA_OFFSET_L, data[3]);
  653.   writeByte(MPU9250_ADDRESS, ZA_OFFSET_H, data[4]);
  654.   writeByte(MPU9250_ADDRESS, ZA_OFFSET_L, data[5]);

  656.   // Output scaled accelerometer biases for display in the main program
  657.   dest2[0] = (float)accel_bias[0] / (float)accelsensitivity;
  658.   dest2[1] = (float)accel_bias[1] / (float)accelsensitivity;
  659.   dest2[2] = (float)accel_bias[2] / (float)accelsensitivity;
  660. }

  662. // Accelerometer and gyroscope self test; check calibration wrt factory settings
  663. void MPU9250SelfTest(float * destination) // Should return percent deviation from factory trim values, +/- 14 or less deviation is a pass
  664. {
  665.   uint8_t rawData[6] = {0, 0, 0, 0, 0, 0};
  666.   uint8_t selfTest[6];
  667.   int32_t gAvg[3] = {0}, aAvg[3] = {0}, aSTAvg[3] = {0}, gSTAvg[3] = {0};
  668.   float factoryTrim[6];
  669.   uint8_t FS = 0;

  671.   writeByte(MPU9250_ADDRESS, SMPLRT_DIV, 0x00); // Set gyro sample rate to 1 kHz
  672.   writeByte(MPU9250_ADDRESS, CONFIG, 0x02); // Set gyro sample rate to 1 kHz and DLPF to 92 Hz
  673.   writeByte(MPU9250_ADDRESS, GYRO_CONFIG, FS << 3); // Set full scale range for the gyro to 250 dps
  674.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG2, 0x02); // Set accelerometer rate to 1 kHz and bandwidth to 92 Hz
  675.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG, FS << 3); // Set full scale range for the accelerometer to 2 g

  677.   for ( int ii = 0; ii < 200; ii++) { // get average current values of gyro and acclerometer

  679.     readBytes(MPU9250_ADDRESS, ACCEL_XOUT_H, 6, &rawData[0]); // Read the six raw data registers into data array
  680.     aAvg[0] += (int16_t)(((int16_t)rawData[0] << 8) | rawData[1]) ; // Turn the MSB and LSB into a signed 16-bit value
  681.     aAvg[1] += (int16_t)(((int16_t)rawData[2] << 8) | rawData[3]) ;
  682.     aAvg[2] += (int16_t)(((int16_t)rawData[4] << 8) | rawData[5]) ;

  684.     readBytes(MPU9250_ADDRESS, GYRO_XOUT_H, 6, &rawData[0]); // Read the six raw data registers sequentially into data array
  685.     gAvg[0] += (int16_t)(((int16_t)rawData[0] << 8) | rawData[1]) ; // Turn the MSB and LSB into a signed 16-bit value
  686.     gAvg[1] += (int16_t)(((int16_t)rawData[2] << 8) | rawData[3]) ;
  687.     gAvg[2] += (int16_t)(((int16_t)rawData[4] << 8) | rawData[5]) ;
  688.   }

  690.   for (int ii = 0; ii < 3; ii++) { // Get average of 200 values and store as average current readings
  691.     aAvg[ii] /= 200;
  692.     gAvg[ii] /= 200;
  693.   }

  695.   // Configure the accelerometer for self-test
  696.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG, 0xE0); // Enable self test on all three axes and set accelerometer range to +/- 2 g
  697.   writeByte(MPU9250_ADDRESS, GYRO_CONFIG, 0xE0); // Enable self test on all three axes and set gyro range to +/- 250 degrees/s
  698.   delay(25); // Delay a while to let the device stabilize

  700.   for ( int ii = 0; ii < 200; ii++) { // get average self-test values of gyro and acclerometer

  702.     readBytes(MPU9250_ADDRESS, ACCEL_XOUT_H, 6, &rawData[0]); // Read the six raw data registers into data array
  703.     aSTAvg[0] += (int16_t)(((int16_t)rawData[0] << 8) | rawData[1]) ; // Turn the MSB and LSB into a signed 16-bit value
  704.     aSTAvg[1] += (int16_t)(((int16_t)rawData[2] << 8) | rawData[3]) ;
  705.     aSTAvg[2] += (int16_t)(((int16_t)rawData[4] << 8) | rawData[5]) ;

  707.     readBytes(MPU9250_ADDRESS, GYRO_XOUT_H, 6, &rawData[0]); // Read the six raw data registers sequentially into data array
  708.     gSTAvg[0] += (int16_t)(((int16_t)rawData[0] << 8) | rawData[1]) ; // Turn the MSB and LSB into a signed 16-bit value
  709.     gSTAvg[1] += (int16_t)(((int16_t)rawData[2] << 8) | rawData[3]) ;
  710.     gSTAvg[2] += (int16_t)(((int16_t)rawData[4] << 8) | rawData[5]) ;
  711.   }

  713.   for (int ii = 0; ii < 3; ii++) { // Get average of 200 values and store as average self-test readings
  714.     aSTAvg[ii] /= 200;
  715.     gSTAvg[ii] /= 200;
  716.   }

  718.   // Configure the gyro and accelerometer for normal operation
  719.   writeByte(MPU9250_ADDRESS, ACCEL_CONFIG, 0x00);
  720.   writeByte(MPU9250_ADDRESS, GYRO_CONFIG, 0x00);
  721.   delay(25); // Delay a while to let the device stabilize

  723.   // Retrieve accelerometer and gyro factory Self-Test Code from USR_Reg
  724.   selfTest[0] = readByte(MPU9250_ADDRESS, SELF_TEST_X_ACCEL); // X-axis accel self-test results
  725.   selfTest[1] = readByte(MPU9250_ADDRESS, SELF_TEST_Y_ACCEL); // Y-axis accel self-test results
  726.   selfTest[2] = readByte(MPU9250_ADDRESS, SELF_TEST_Z_ACCEL); // Z-axis accel self-test results
  727.   selfTest[3] = readByte(MPU9250_ADDRESS, SELF_TEST_X_GYRO); // X-axis gyro self-test results
  728.   selfTest[4] = readByte(MPU9250_ADDRESS, SELF_TEST_Y_GYRO); // Y-axis gyro self-test results
  729.   selfTest[5] = readByte(MPU9250_ADDRESS, SELF_TEST_Z_GYRO); // Z-axis gyro self-test results

  731.   // Retrieve factory self-test value from self-test code reads
  732.   factoryTrim[0] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[0] - 1.0) )); // FT[Xa] factory trim calculation
  733.   factoryTrim[1] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[1] - 1.0) )); // FT[Ya] factory trim calculation
  734.   factoryTrim[2] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[2] - 1.0) )); // FT[Za] factory trim calculation
  735.   factoryTrim[3] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[3] - 1.0) )); // FT[Xg] factory trim calculation
  736.   factoryTrim[4] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[4] - 1.0) )); // FT[Yg] factory trim calculation
  737.   factoryTrim[5] = (float)(2620 / 1 << FS) * (pow( 1.01 , ((float)selfTest[5] - 1.0) )); // FT[Zg] factory trim calculation

  739.   // Report results as a ratio of (STR - FT)/FT; the change from Factory Trim of the Self-Test Response
  740.   // To get percent, must multiply by 100
  741.   for (int i = 0; i < 3; i++) {
  742.     destination[i] = 100.0 * ((float)(aSTAvg[i] - aAvg[i])) / factoryTrim[i] - 100.; // Report percent differences
  743.     destination[i + 3] = 100.0 * ((float)(gSTAvg[i] - gAvg[i])) / factoryTrim[i + 3] - 100.; // Report percent differences
  744.   }

  746. }



  750. void writeByte(uint8_t address, uint8_t subAddress, uint8_t data)
  751. {
  752.   Wire1.beginTransmission(address); // Initialize the Tx buffer
  753.   Wire1.write(subAddress); // Put slave register address in Tx buffer
  754.   Wire1.write(data); // Put data in Tx buffer
  755.   Wire1.endTransmission(); // Send the Tx buffer
  756. }

  758. uint8_t readByte(uint8_t address, uint8_t subAddress)
  759. {
  760.   uint8_t data; // `data` will store the register data
  761.   Wire1.beginTransmission(address); // Initialize the Tx buffer
  762.   Wire1.write(subAddress);     // Put slave register address in Tx buffer
  763.   Wire1.endTransmission(I2C_NOSTOP); // Send the Tx buffer, but send a restart to keep connection alive
  764.   //    Wire1.endTransmission(false); // Send the Tx buffer, but send a restart to keep connection alive
  765.   //    Wire1.requestFrom(address, 1); // Read one byte from slave register address
  766.   Wire1.requestFrom(address, (size_t) 1); // Read one byte from slave register address
  767.   data = Wire1.read(); // Fill Rx buffer with result
  768.   return data; // Return data read from slave register
  769. }

  771. void readBytes(uint8_t address, uint8_t subAddress, uint8_t count, uint8_t * dest)
  772. {
  773.   Wire1.beginTransmission(address); // Initialize the Tx buffer
  774.   Wire1.write(subAddress); // Put slave register address in Tx buffer
  775.   Wire1.endTransmission(I2C_NOSTOP); // Send the Tx buffer, but send a restart to keep connection alive
  776.   //    Wire1.endTransmission(false); // Send the Tx buffer, but send a restart to keep connection alive
  777.   uint8_t i = 0;
  778.   // Wire1.requestFrom(address, count); // Read bytes from slave register address
  779.   Wire1.requestFrom(address, (size_t) count); // Read bytes from slave register address
  780.   while (Wire1.available()) {
  781.     dest[i++] = Wire1.read();
  782.   } // Put read results in the Rx buffer
  783. }


点击(此处)折叠或打开

  1. // Implementation of Sebastian Madgwick's "...efficient orientation filter for... inertial/magnetic sensor arrays"
  2. // (see http://www.x-io.co.uk/category/open-source/ for examples and more details)
  3. // which fuses acceleration, rotation rate, and magnetic moments to produce a quaternion-based estimate of absolute
  4. // device orientation -- which can be converted to yaw, pitch, and roll. Useful for stabilizing quadcopters, etc.
  5. // The performance of the orientation filter is at least as good as conventional Kalman-based filtering algorithms
  6. // but is much less computationally intensive---it can be performed on a 3.3 V Pro Mini operating at 8 MHz!
  7.         void MadgwickQuaternionUpdate(float ax, float ay, float az, float gx, float gy, float gz, float mx, float my, float mz)
  8.         {
  9.             float q1 = q[0], q2 = q[1], q3 = q[2], q4 = q[3]; // short name local variable for readability
  10.             float norm;
  11.             float hx, hy, _2bx, _2bz;
  12.             float s1, s2, s3, s4;
  13.             float qDot1, qDot2, qDot3, qDot4;

  15.             // Auxiliary variables to avoid repeated arithmetic
  16.             float _2q1mx;
  17.             float _2q1my;
  18.             float _2q1mz;
  19.             float _2q2mx;
  20.             float _4bx;
  21.             float _4bz;
  22.             float _2q1 = 2.0f * q1;
  23.             float _2q2 = 2.0f * q2;
  24.             float _2q3 = 2.0f * q3;
  25.             float _2q4 = 2.0f * q4;
  26.             float _2q1q3 = 2.0f * q1 * q3;
  27.             float _2q3q4 = 2.0f * q3 * q4;
  28.             float q1q1 = q1 * q1;
  29.             float q1q2 = q1 * q2;
  30.             float q1q3 = q1 * q3;
  31.             float q1q4 = q1 * q4;
  32.             float q2q2 = q2 * q2;
  33.             float q2q3 = q2 * q3;
  34.             float q2q4 = q2 * q4;
  35.             float q3q3 = q3 * q3;
  36.             float q3q4 = q3 * q4;
  37.             float q4q4 = q4 * q4;

  39.             // Normalise accelerometer measurement
  40.             norm = sqrtf(ax * ax + ay * ay + az * az);
  41.             if (norm == 0.0f) return; // handle NaN
  42.             norm = 1.0f/norm;
  43.             ax *= norm;
  44.             ay *= norm;
  45.             az *= norm;

  47.             // Normalise magnetometer measurement
  48.             norm = sqrtf(mx * mx + my * my + mz * mz);
  49.             if (norm == 0.0f) return; // handle NaN
  50.             norm = 1.0f/norm;
  51.             mx *= norm;
  52.             my *= norm;
  53.             mz *= norm;

  55.             // Reference direction of Earth's magnetic field
  56.             _2q1mx = 2.0f * q1 * mx;
  57.             _2q1my = 2.0f * q1 * my;
  58.             _2q1mz = 2.0f * q1 * mz;
  59.             _2q2mx = 2.0f * q2 * mx;
  60.             hx = mx * q1q1 - _2q1my * q4 + _2q1mz * q3 + mx * q2q2 + _2q2 * my * q3 + _2q2 * mz * q4 - mx * q3q3 - mx * q4q4;
  61.             hy = _2q1mx * q4 + my * q1q1 - _2q1mz * q2 + _2q2mx * q3 - my * q2q2 + my * q3q3 + _2q3 * mz * q4 - my * q4q4;
  62.             _2bx = sqrtf(hx * hx + hy * hy);
  63.             _2bz = -_2q1mx * q3 + _2q1my * q2 + mz * q1q1 + _2q2mx * q4 - mz * q2q2 + _2q3 * my * q4 - mz * q3q3 + mz * q4q4;
  64.             _4bx = 2.0f * _2bx;
  65.             _4bz = 2.0f * _2bz;

  67.             // Gradient decent algorithm corrective step
  68.             s1 = -_2q3 * (2.0f * q2q4 - _2q1q3 - ax) + _2q2 * (2.0f * q1q2 + _2q3q4 - ay) - _2bz * q3 * (_2bx * (0.5f - q3q3 - q4q4) + _2bz * (q2q4 - q1q3) - mx) + (-_2bx * q4 + _2bz * q2) * (_2bx * (q2q3 - q1q4) + _2bz * (q1q2 + q3q4) - my) + _2bx * q3 * (_2bx * (q1q3 + q2q4) + _2bz * (0.5f - q2q2 - q3q3) - mz);
  69.             s2 = _2q4 * (2.0f * q2q4 - _2q1q3 - ax) + _2q1 * (2.0f * q1q2 + _2q3q4 - ay) - 4.0f * q2 * (1.0f - 2.0f * q2q2 - 2.0f * q3q3 - az) + _2bz * q4 * (_2bx * (0.5f - q3q3 - q4q4) + _2bz * (q2q4 - q1q3) - mx) + (_2bx * q3 + _2bz * q1) * (_2bx * (q2q3 - q1q4) + _2bz * (q1q2 + q3q4) - my) + (_2bx * q4 - _4bz * q2) * (_2bx * (q1q3 + q2q4) + _2bz * (0.5f - q2q2 - q3q3) - mz);
  70.             s3 = -_2q1 * (2.0f * q2q4 - _2q1q3 - ax) + _2q4 * (2.0f * q1q2 + _2q3q4 - ay) - 4.0f * q3 * (1.0f - 2.0f * q2q2 - 2.0f * q3q3 - az) + (-_4bx * q3 - _2bz * q1) * (_2bx * (0.5f - q3q3 - q4q4) + _2bz * (q2q4 - q1q3) - mx) + (_2bx * q2 + _2bz * q4) * (_2bx * (q2q3 - q1q4) + _2bz * (q1q2 + q3q4) - my) + (_2bx * q1 - _4bz * q3) * (_2bx * (q1q3 + q2q4) + _2bz * (0.5f - q2q2 - q3q3) - mz);
  71.             s4 = _2q2 * (2.0f * q2q4 - _2q1q3 - ax) + _2q3 * (2.0f * q1q2 + _2q3q4 - ay) + (-_4bx * q4 + _2bz * q2) * (_2bx * (0.5f - q3q3 - q4q4) + _2bz * (q2q4 - q1q3) - mx) + (-_2bx * q1 + _2bz * q3) * (_2bx * (q2q3 - q1q4) + _2bz * (q1q2 + q3q4) - my) + _2bx * q2 * (_2bx * (q1q3 + q2q4) + _2bz * (0.5f - q2q2 - q3q3) - mz);
  72.             norm = sqrtf(s1 * s1 + s2 * s2 + s3 * s3 + s4 * s4); // normalise step magnitude
  73.             norm = 1.0f/norm;
  74.             s1 *= norm;
  75.             s2 *= norm;
  76.             s3 *= norm;
  77.             s4 *= norm;

  79.             // Compute rate of change of quaternion
  80.             qDot1 = 0.5f * (-q2 * gx - q3 * gy - q4 * gz) - beta * s1;
  81.             qDot2 = 0.5f * (q1 * gx + q3 * gz - q4 * gy) - beta * s2;
  82.             qDot3 = 0.5f * (q1 * gy - q2 * gz + q4 * gx) - beta * s3;
  83.             qDot4 = 0.5f * (q1 * gz + q2 * gy - q3 * gx) - beta * s4;

  85.             // Integrate to yield quaternion
  86.             q1 += qDot1 * deltat;
  87.             q2 += qDot2 * deltat;
  88.             q3 += qDot3 * deltat;
  89.             q4 += qDot4 * deltat;
  90.             norm = sqrtf(q1 * q1 + q2 * q2 + q3 * q3 + q4 * q4); // normalise quaternion
  91.             norm = 1.0f/norm;
  92.             q[0] = q1 * norm;
  93.             q[1] = q2 * norm;
  94.             q[2] = q3 * norm;
  95.             q[3] = q4 * norm;

  97.         }
  98.   
  99.   
  100.   
  101.  // Similar to Madgwick scheme but uses proportional and integral filtering on the error between estimated reference vectors and
  102.  // measured ones.
  103.             void MahonyQuaternionUpdate(float ax, float ay, float az, float gx, float gy, float gz, float mx, float my, float mz)
  104.         {
  105.             float q1 = q[0], q2 = q[1], q3 = q[2], q4 = q[3]; // short name local variable for readability
  106.             float norm;
  107.             float hx, hy, bx, bz;
  108.             float vx, vy, vz, wx, wy, wz;
  109.             float ex, ey, ez;
  110.             float pa, pb, pc;

  112.             // Auxiliary variables to avoid repeated arithmetic
  113.             float q1q1 = q1 * q1;
  114.             float q1q2 = q1 * q2;
  115.             float q1q3 = q1 * q3;
  116.             float q1q4 = q1 * q4;
  117.             float q2q2 = q2 * q2;
  118.             float q2q3 = q2 * q3;
  119.             float q2q4 = q2 * q4;
  120.             float q3q3 = q3 * q3;
  121.             float q3q4 = q3 * q4;
  122.             float q4q4 = q4 * q4;

  124.             // Normalise accelerometer measurement
  125.             norm = sqrtf(ax * ax + ay * ay + az * az);
  126.             if (norm == 0.0f) return; // handle NaN
  127.             norm = 1.0f / norm; // use reciprocal for division
  128.             ax *= norm;
  129.             ay *= norm;
  130.             az *= norm;

  132.             // Normalise magnetometer measurement
  133.             norm = sqrtf(mx * mx + my * my + mz * mz);
  134.             if (norm == 0.0f) return; // handle NaN
  135.             norm = 1.0f / norm; // use reciprocal for division
  136.             mx *= norm;
  137.             my *= norm;
  138.             mz *= norm;

  140.             // Reference direction of Earth




采用 MPU9250的话.
首先配置
MPU_9250_I2C_SLAVE_ADDR = 0x68        MPU_9250_I2C_DEVICE_ID = 0x71

  1.         write_byte_data(PWR_MGMT_1, 0x01)  /* Reset device. */
  2.         time_sleep(0.1)
  3.         write_byte_data(PWR_MGMT_1, 0x00)  /* Wake up chip. */

  5.         write_byte_data(CONFIG, 0x03)       //0x1A --> accel bandwidth=44 delay=4.9ms  gyroscope:band width=42Hz, delay=4.8ms 1Khz  : 低通滤波
  6.         write_byte_data(SMPLRT_DIV, 0x09)    //0x19  --> 采样频率= 陀螺仪输出频率 / （ 1+SMPLRT_DIV).  根据低通滤波设置, 陀螺仪输出频率为 1KHz. 这样设定采样频率为 100Hz. 即 FIFO输出/DMP采样/运动检测都基于此采样频率.
  7.         write_byte_data(GYRO_CONFIG, gyro_range << 3) //0x1B--> GYRO_RANGE_250
  8.         write_byte_data(ACCEL_CONFIG, accel_range << 3) //0x1C --> ACCEL_RANGE_2G
  9.         write_byte_data(ACCEL_CONFIG_2, 0x03) //低通滤波二配置.
  10.         write_byte_data(INT_PIN_CFG, 0x02) //0x37 --> 直接访问扩展IIC.即访问磁力计
从寄存器的取值为 ADC的值.  16 位 ADC 的 加速度/ 角速度 对应的实际数值转化过来, 扣除最初的符号位. 
实际加速度值 = 寄存器值 / (精度 2.0g/32768.0)  // ACCEL_RANGE_2G 时
实际角速度值 = (寄存器值 / (精度 250.0/32768.0) )  * PI/180.0f     // GYRO_RANGE_250 时
最高为1时, 记得变换成 负值.

角速度

采样完毕, INT_STATUS = 0x3A 寄存器 0x01 设置, 从寄存器 ACCEL_OUT 0x3B/GYRO_OUT 0x41 读出信息, 转化成实际速度值.

计算四元数


点击(此处)折叠或打开

  1. #include "IMU.h"

  3. #include "math.h"

  5. #include "Maths.h"

  7. // ==================================================================================

  9. //加Q： 4 4 8 2 2 2 46 5 讨论

  11. // 描述:

  13. // 必须定义'halfT '为周期的一半，以及滤波器的参数Kp和Ki

  15. // 四元数'q0', 'q1', 'q2', 'q3'定义为全局变量

  17. // 需要在每一个采样周期调用'IMUupdate()'函数

  19. // 陀螺仪数据单位是弧度/秒，加速度计的单位无关重要，因为会被规范化

  21. // ==================================================================================

  23. #define Kp 90.0f // 比例常数

  25. #define Ki 0.001f // 积分常数

  27. #define halfT 0.0005f//半周期

  29. #define T    0.001f // 周期为1ms

  31. // ==================================================================================

  33. // 变量定义

  35. float q0 = 1, q1 = 0, q2 = 0, q3 = 0; // 四元数

  37. float exInt = 0, eyInt = 0, ezInt = 0; // 误差积分累计值

  39. // ==================================================================================

  41. // 函数原型：void IMUupdate(float gx, float gy, float gz, float ax, float ay, float az)

  43. // 功 能：互补滤波进行姿态解算

  45. // 输 入：陀螺仪数据及加速度计数据

  47. // ==================================================================================

  49. void IMUupdate(float gx, float gy, float gz, float ax, float ay, float az)

  51. {

  53. float norm;

  55. float vx, vy, vz;

  57. float ex, ey, ez;

  59. //四元数积分，求得当前的姿态

  61. float q0_last = q0;

  63. float q1_last = q1;

  65. float q2_last = q2;

  67. float q3_last = q3;

  69. //把加速度计的三维向量转成单位向量

  71. norm = invSqrt(ax*ax + ay*ay + az*az);

  73. ax = ax * norm;

  75. ay = ay * norm;

  77. az = az * norm;

  79. //估计重力加速度方向在飞行器坐标系中的表示，为四元数表示的旋转矩阵的第三行

  81. vx = 2*(q1*q3 - q0*q2);

  83. vy = 2*(q0*q1 + q2*q3);

  85. vz = q0*q0 - q1*q1 - q2*q2 + q3*q3;

  87. //加速度计读取的方向与重力加速度方向的差值，用向量叉乘计算

  89. ex = ay*vz - az*vy;

  91. ey = az*vx - ax*vz;

  93. ez = ax*vy - ay*vx;

  95. //误差累积，已与积分常数相乘

  97. exInt = exInt + ex*Ki;

  99. eyInt = eyInt + ey*Ki;

  101. ezInt = ezInt + ez*Ki;

  103. //用叉积误差来做PI修正陀螺零偏，即抵消陀螺读数中的偏移量

  105. gx = gx + Kp*ex + exInt;

  107. gy = gy + Kp*ey + eyInt;

  109. gz = gz + Kp*ez + ezInt;

  111. //一阶近似算法

  113. q0 = q0_last + (-q1_last*gx - q2_last*gy - q3_last*gz)*halfT;

  115. q1 = q1_last + ( q0_last*gx + q2_last*gz - q3_last*gy)*halfT;

  117. q2 = q2_last + ( q0_last*gy - q1_last*gz + q3_last*gx)*halfT;

  119. q3 = q3_last + ( q0_last*gz + q1_last*gy - q2_last*gx)*halfT;

  121. //    //二阶近似算法

  123. //    float delta2 = (gx*gx + gy*gy + gz*gz)*T*T;

  125. //    q0 = q0_last*(1-delta2/8) + (-q1_last*gx - q2_last*gy - q3_last*gz)*halfT;

  127. //    q1 = q1_last*(1-delta2/8) + ( q0_last*gx + q2_last*gz - q3_last*gy)*halfT;

  129. //    q2 = q2_last*(1-delta2/8) + ( q0_last*gy - q1_last*gz + q3_last*gx)*halfT;

  131. //    q3 = q3_last*(1-delta2/8) + ( q0_last*gz + q1_last*gy - q2_last*gx)*halfT;

  133. //    //三阶近似算法

  135. //    float delta2 = (gx*gx + gy*gy + gz*gz)*T*T;

  137. //    q0 = q0_last*(1-delta2/8) + (-q1_last*gx - q2_last*gy - q3_last*gz)*T*(0.5 - delta2/48);

  139. //    q1 = q1_last*(1-delta2/8) + ( q0_last*gx + q2_last*gz - q3_last*gy)*T*(0.5 - delta2/48);

  141. //    q2 = q2_last*(1-delta2/8) + ( q0_last*gy - q1_last*gz + q3_last*gx)*T*(0.5 - delta2/48);

  143. //    q3 = q3_last*(1-delta2/8) + ( q0_last*gz + q1_last*gy - q2_last*gx)*T*(0.5 - delta2/48);

  145. //    //四阶近似算法

  147. //    float delta2 = (gx*gx + gy*gy + gz*gz)*T*T;

  149. //    q0 = q0_last*(1 - delta2/8 + delta2*delta2/384) + (-q1_last*gx - q2_last*gy - q3_last*gz)*T*(0.5 - delta2/48);

  151. //    q1 = q1_last*(1 - delta2/8 + delta2*delta2/384) + ( q0_last*gx + q2_last*gz - q3_last*gy)*T*(0.5 - delta2/48);

  153. //    q2 = q2_last*(1 - delta2/8 + delta2*delta2/384) + ( q0_last*gy - q1_last*gz + q3_last*gx)*T*(0.5 - delta2/48);

  155. //    q3 = q3_last*(1 - delta2/8 + delta2*delta2/384) + ( q0_last*gz + q1_last*gy - q2_last*gx)*T*(0.5 - delta2/48);

  157. //四元数规范化

  159. norm = invSqrt(q0*q0 + q1*q1 + q2*q2 + q3*q3);

  161. q0 = q0 * norm;

  163. q1 = q1 * norm;

  165. q2 = q2 * norm;

  167. q3 = q3 * norm;

  169. //out_angle.yaw += filter_gyro.z * RawData_to_Angle * 0.001f;

  171. out_angle.yaw = atan2(2 * q1 * q2 + 2 * q0 * q3, -2 * q2*q2 - 2 * q3* q3 + 1)* 57.3;

  173. }

  175. //Pitch是围绕X轴旋转，也叫作俯仰角

  177. //roll是围绕Y轴旋转，也叫作翻滚角

  179. //yaw是围绕Z轴旋转，也叫作偏航角

  181. /******************************************************************************

  183. 函数原型：    void Get_Eulerian_Angle(struct _out_angle *angle)

  185. 功 能：    四元数转欧拉角

  187. *******************************************************************************/

  189. void Get_Eulerian_Angle(struct _out_angle *angle)

  191. {

  193. angle->pitch = -atan2(2.0f*(q0*q1 + q2*q3),q0*q0 - q1*q1 - q2*q2 + q3*q3)*Radian_to_Angle;

  195. angle->roll = asin (2.0f*(q0*q2 - q1*q3))*Radian_to_Angle;

  197. Angle.pitch = angle->pitch;

  199. Angle.roll = angle->roll;

  201. Angle.yaw = out_angle.yaw;

  203. }




阅读(10265) | 评论(0) | 转发(0) |
0

上一篇：MDK KEIL中使用STM32F4XX芯片硬件浮点单元FPU

下一篇：用加速度和角速度计算位移

给主人留下些什么吧！~~
关于我们 | 关于IT168 | 联系方式 | 广告合作 | 法律声明 | 免费注册

Copyright 2001-2010 ChinaUnix.net All Rights Reserved 北京皓辰网域网络信息技术有限公司. 版权所有

感谢所有关心和支持过ChinaUnix的朋友们

16024965号-6