Research Article  Open Access
Indoor Dead Reckoning Localization Using Ultrasonic Anemometer with IMU
Abstract
Dead reckoning is an important aspect of estimating the instantaneous position of a mobile robot. An inertial measurement unit (IMU) is generally used for dead reckoning because it measures triaxis acceleration and triaxis angular velocities in order to estimate the position of the mobile robot. Positioning with inertial data is reasonable for a short period of time. However, the velocity, position, and attitude errors increase over time. Much research has been conducted in ways to reduce these errors. To position a mobile robot, an absolute positioning method can be combined with dead reckoning. The performance of a combined positioning method can be improved based on improvement in dead reckoning. In this paper, an ultrasonic anemometer is used to improve the performance of dead reckoning when indoors. A new approach to the equation of an ultrasonic anemometer is proposed. The ultrasonic anemometer prevents divergence of the mobile robot’s velocity. To position a mobile robot indoors, the ultrasonic anemometer measures the relative movement of air while the robot moves through static air. Velocity data from the ultrasonic anemometer and the acceleration and angular velocity data from the IMU are combined via Kalman filter. Finally we show that the proposed method has the performance with a positioning method using encoders on a good floor condition.
1. Introduction
Study on localization, which is an important technique in order to achieve the purposes of a mobile robot, has become popular [1–20]. There are two methods of localization [1–3]. One is the absolute positioning method, which uses external devices that include absolute position information or a direction measuring system. The other is the relative positioning method, which mostly uses inertial sensors or encoders. Generally, two or more methods are combined to more accurately estimate position [14–17].
The absolute positioning method for outdoor applications uses a global navigation satellite system (GNSS). The GNSS receiver calculates position using data received from several GNSS satellite. This method has a low update rate and cannot be used for indoor system because the signal may be too weak or not reachable to detect the signal [2, 14]. For indoor applications, an ultrasonic local positioning system, an infrared network system, or radio frequency identification (RFID) is commonly used. The ultrasonic local positioning system, infrared network system, and RFID usually measure distances or differences in distance by time of flight (TOF), time difference of flight (TDOF), or received signal strength (RSS). Some of these systems measure angle of arrival (AOA) from transmitter to receiver. These systems cannot measure long distances, require additional equipment, and are affected by signal interference. The absolute positioning methods do not accumulate positioning errors, but the positioning errors could be relatively large. Sometimes, a laser range finder (LRF) is used for localization of a mobile robot [20]. The positioning method using LRF uses a simultaneous localization and mapbuilding (SLAM) method or triangulation method. The SLAM method measures the distance and direction of a wall or an obstacle, creates a map, and calculates the position of the mobile robot from the map. The triangulation method is used with a camera or a reflector to calculate the position of the mobile robot by measuring the distance or direction of the position where the position is known in advance. An LRF can measure long distances, provide highaccuracy range data, and does not need external equipment. But it is used for very expensive mobile robots because of its price.
The relative positioning method uses an inertial measurement unit (IMU), which consists of a triaxis accelerometer and a triaxis gyroscope, or encoders attached to the mobile robot’s wheels [1–7]. The inertial data from the IMU are accumulated in order to estimate position, velocity, and attitude. The velocity data of each wheel are estimated by counting pulses from each encoder. These methods are simple, lowcost, and easy to implement. Each error in a relative positioning system is small, but these errors accumulate during calculations of position, velocity, and attitude.
A lot of researches into positioning have used the IMU [1–13]. There are accumulated errors during integrative calculations to find velocity or position. Reducing error accumulation is an important technique when using an IMU. If the performance of dead reckoning is improved, the performance of an integrated positioning system with an absolute positioning method can also be improved. Some researchers try to reduce accumulated localization or attitude errors by using IMU data [8–13]. For localization errors, some researchers use zero velocity compensation to remove accumulated errors [8]. The zero velocity compensation method cannot be applied to a continuously moving mobile robot. Other researchers adopt odometry to compensate for accumulated velocity errors [1–7], but the errors from odometry can be large in some situations [2, 3]; for example, the odometry error becomes large when a mobile robot moves on a slippery or uneven floor.
In this paper, an ultrasonic anemometer is adopted in order to measure the velocity of a mobile robot. Generally, the ultrasonic anemometer measures wind velocity and direction at a fixed location. But, an ultrasonic anemometer can measure relative velocity of air while the mobile robot moves through static air. A new approach to the equation of an ultrasonic anemometer is proposed. The traditional ultrasonic anemometer equations contain logical errors and derive new equations considering this. Data from the ultrasonic anemometer will include noise, because there is turbulent air flow, and the mobile robot uses only a small measuring range of the ultrasonic anemometer. But these errors can be compensated for by combining data from odometry or an IMU.
The ultrasonic anemometer can be used without stationary motion and can measure the velocity of a mobile robot directly as a method of localization from dead reckoning. Also a Kalman filter is designed to combine data from the accelerometer and the ultrasonic anemometer. The proposed system provides a technique for designing an autonomous mobile robot at a reasonable price. The performance of the proposed system was evaluated with experiments.
The rest of this paper is organized as follows. Section 2 describes the theoretical background of each localization method. Section 3 describes the proposed localization algorithm using an ultrasonic anemometer and an IMU. The experimental results are shown in Section 4. The final conclusions are in Section 5.
2. The Navigation System
2.1. The Ultrasonic Anemometer
The proposed system uses a continuouswave ultrasonic anemometer to measure the velocity of a mobile robot. The ultrasonic anemometer measures relative velocity of air while the mobile robot moves through static air.
Figure 1 shows a simplified model of an ultrasonic anemometer. The system is simplified as one transmitter and one receiver to solve the relationship between the velocity of the mobile robot and the flight time of ultrasound. The traditional equation of an ultrasonic anemometer is derived from the following equation:
If the angle of the ultrasonic wave is perpendicular to the moving direction of the air, the following is obtained:
The above equation has a constant value regardless of wind speed. The movement of air when an ultrasonic anemometer is stationary can be seen to be the same as the movement of an ultrasonic anemometer in stationary air. Therefore, it is reasonable to derive the equation as follows:
The mobile robot moves with velocity . The transmitter generates ultrasound waves; the mobile robot moves distance , and the receiver receives the ultrasound waves, which traverse distance at the velocity of sound, . The flight time of ultrasound is calculated as follows.
The simplified model can measure the magnitude of the velocity of a mobile robot. The ultrasonic anemometer has to have two receivers to measure the magnitude and direction of velocity in one dimension, as seen in Figure 2.
The flight time of ultrasonic waves from transmitter to receiver 1, , and receiver 2, , is calculated as follows:
The ultrasonic anemometer uses time difference of flight for two receivers
The traditional equation corresponding to (5) is as follows:
The results show that the proposed equation is simpler than the traditional equation.
The velocity of mobile robot can be calculated from time difference as follows:
The approximation can be applied when the speed of the mobile robot is much smaller than the speed of sound.
The processor measures the phase difference from each receiver, whereas the transmitter generates continuous ultrasound waves.
2.2. The Inertial Measurement Unit
The IMU consists of a triaxial accelerometer and a triaxial gyroscope and sometimes also has a triaxial or twinaxial magnetometer. The triaxial accelerometer measures threedimension acceleration. The inclination of the mobile robot is estimated from the triaxial accelerometer when the mobile robot is stationary. The triaxial gyroscope measures the threedimension angular rate of the mobile robot. The magnetometer is used for accumulated attitude error compensation in the mobile robot.
A calibration procedure is needed before using the inertial sensors, which are the accelerometer and the gyroscope. For an inertial navigation system, the positon of the mobile robot is estimated by accumulating inertial data from the initial position and the initial attitude. Therefore, sensor bias and the initial attitude can cause serious position errors. The origin is assumed to be the initial position. Sensor calibration is used to calculate the initial attitude and cancel out the bias errors of the inertial sensors when the mobile robot is stationary. The sensor calibration process is depicted in Figure 3.
When the mobile robot is stationary the acceleration on body frame is as follows:
is the acceleration on the navigation frame, is the gravitational acceleration, is the pitch angle, and is the roll angle.
is the transform matrix from the navigation frame to the body frame
The notation is a cosine and the notation is a sine. is the yaw angle of the mobile robot. The twinaxial magnetometer measures the yaw angle which is based on magnetic north.
The roll angle and pitch angle, which are from (8), are as follows:
The measurement of the accelerometer is
is the transpose matrix of .
When the mobile robot is stationary, the bias of the accelerometer, , and the estimated acceleration in the body frame, , are
The gyroscope data also have to be calibrated with the following equations:
A Kalman filter is adopted to estimate the position of the mobile robot, which is described in the next section.
3. The Fusion of Ultrasonic Anemometer and IMU
The Kalman filter is designed to combine data from the IMU, ultrasonic anemometer, and magnetometer. The Kalman filter is a kind of optimal stochastic filter, which has reasonable performance and low complexity and is suitable for realtime systems. The proposed system is modeled as follows:
State variable is a linear combination of previous state variable and input . Measurement is a linear combination of state variable . For both equations, the process noise, , and the measurement noise, , are added.
The Kalman filter consists of a prediction step and a correction step. The former predicts the a priori state estimation from the previous a posteriori state estimation and input data. The latter updates measurements with the a posteriori state estimation and measurement data. A block diagram of the Kalman filter is depicted in Figure 4.
The state variable is defined for the Kalman filter as follows:
Variable is the velocity of the mobile robot in the navigation frame, and variable is the attitude of the mobile robot in the navigation frame.
The predictor equations for the Kalman filter are as follows:
The a priori estimation of the state variable is , and the estimation error covariance matrix is . The error covariance matrix of process error is .
Variable is the prediction input from the accelerometers and the gyroscope. Matrix is the state of the transient matrix, which is an identity matrix. Matrix is the prediction input matrix, which is a rotation matrix from body frame to navigation frame
Matrix is a threebythree zero matrix.
The corrector equations for the Kalman filter are as follows:
The a posteriori estimation of the state variable is , and the estimation error covariance matrix is . The error covariance matrix of correction is . The Kalman gain, , is calculated from the a priori estimation error covariance matrix and correction error covariance matrix . Variable is the correction input from the ultrasonic anemometer and the magnetometer. In order to minimize magnetic interference, this paper applies the Kalman filter correction error covariance matrix not only the output noise of the magnetic sensor but also the magnetic field distortion measured indoor. As a result, if the mobile robot continues to move, the instantaneous attitude becomes a cumulative output of the gyroscope and gradually converges to the output of the magnetometer. Therefore, if the mobile robot continues to move, it will be less affected by the magnetic field distortion, but if the mobile robot is stationary, the attitude of the mobile robot is greatly affected by the distorted magnetic field. Matrix is the correction matrix, which is a rotation matrix from navigation frame to body frame
4. Experimental Results
The proposed system uses an ultrasonic anemometer in order to measure the velocity of the mobile robot. Generally, an ultrasonic anemometer measures wind velocity and direction at a fixed location. But the ultrasonic anemometer can measure the relative velocity of air while a mobile robot moves through static air. The data from the ultrasonic anemometer contains not only white noise but also turbulent indoor air when the ultrasonic anemometer measures the velocity of the mobile robot. The noise characteristics have to be understood before using the ultrasonic anemometer. Noise from the ultrasonic anemometer in the time domain is shown in Figure 5.
The noise is measured for 20 seconds with a 50 Hz sampling rate when the mobile robot is stationary. The two kinds of noise are shown in Figure 5. One is lowfrequency noise because of the turbulent indoor air. The other is white noise. The histogram of the noise has a Gaussianlike shape, shown in Figure 6. The lowfrequency noise affects more than white noise. The power spectral density of the noise from the ultrasonic anemometer is shown in Figure 7. The power of the noise under 1 Hz has over 88% of the total power. A proper algorithm considering the stochastic characteristics of sensor noise can improve localization performance.
Initial attitude and bias calibration is needed for the inertial sensors. The simple minus operation for sensor calibration can cause significant localization errors on an uneven floor. The inertial data are measured for 20 seconds at a 50 Hz sampling rate when the mobile robot is stationary. The upper graphs in Figures 8 and 9 are the raw data from the accelerometer and the gyroscope. The proposed calibration algorithm using the inertial sensors and the magnetometer is applied to the raw inertial data. The lower graphs in Figures 8 and 9 are the calibrated data. The results of sensor calibration show that the sensor biases are calibrated effectively.
The proposed dead reckoning algorithm was evaluated by experiment with a mobile robot. The mobile robot is based on Stella B2 from NTRex Inc., shown in Figure 10. The Stella B2 has two driving DC motors and its own controller. It is controlled by command via universal asynchronous receivertransmitter and has a positioning computer module installed. An inertial measurement unit (ADIS16354, Analog Devices, Inc.) is used to measure angular velocity and acceleration. An HMC6352 magnetometer measures the absolute yaw attitude of the mobile robot. The ultrasonic anemometer with one transmitter and four receivers measures the absolute velocity of the mobile robot. A TMS320F28335 digital signal processor was used to collect data from the sensors and calculate the position of the mobile robot. The Stella B2 controls itself with a motor encoder, and the positioning computer module sends acceleration, deceleration, constant velocity, and stopping and turning commands to the Stella B2 controller. The positioning computer communicates with a PC in order to send realtime positioning data.
The experiments are conducted on the fabric finished floor which is less slippery than other floors. The positioning by encoders is sensitive to experimental environment. The aim of this experiment is that the proposed system shows the performance compatible with odometry in environment where odometry works well. The positioning method using encoder can cause larger errors in practical applications.
The proposed system was evaluated via linear motion and squaremotion experiments. For linear motion, the mobile robot traveled on a straight path with the velocity command in Figure 11(a). The acceleration rate was 1.07 m/s^{2}, deceleration rate was −1.07 m/s^{2}, and constant velocity was 1.07 m/s. Figure 11 shows the estimated velocity of the mobile robot via accelerometer, ultrasonic anemometer, and the proposed system. The experiment for linear motion is conducted 10 times. The proper experimental results which are considered a general situation are selected. The general situation means that the accumulated velocity error of IMU is not too large, small, or close to zero.
(a) Velocity command
(b) Velocity from the encoders
(c) Velocity from the ultrasonic anemometer
(d) Velocity from the accelerometer
(e) Velocity from the proposed system
Figures 11(b)–11(e) show the velocity estimation with each method. The velocity estimation with the encoders in Figure 11(b) shows the most reasonable data. However, this is not velocity of the mobile robot but wheel revolution. The velocity estimation with the ultrasonic anemometer in Figure 11(c) does not have accumulated velocity error. But, errors from turbulent air and white noise are contained. The velocity estimation with the accelerometer in Figure 11(d) has accumulated error. The accumulated error is only affected by sensor noise when the mobile robot is stationary. But, vibration noise from motors and uneven floor is added when the mobile robot moves. This causes large accumulated error. The velocity estimation with the proposed system in Figure 11(e) shows more reasonable data than data from IMU or ultrasonic anemometer. The velocity from the accelerometer has accumulated error. And the velocity from the ultrasonic anemometer contains more noise. The proposed system shows reasonable velocity estimation performance.
The positioning result for linear motion is in Figure 12. The traditional system uses an accelerometer, a gyroscope, and a magnetometer, whereas the proposed system uses the ultrasonic anemometer with them. The traditional system calculates the position by simply integrating the inertial data. The mobile robot moves from and moves to the north with the velocity instruction in Figure 11(a). The blue line is the trajectory of the position estimation from each method. The actual last position of the mobile robot is marked with a red point.
(a) Positioning by the traditional system
(b) Positioning by the encoders
(c) Positioning by the proposed system
The result of the traditional method in Figure 12(a) shows a longer distance because the error accumulation of an accelerometer is positive. It is hard to estimate the position of the mobile robot using only IMU because the mechanical vibration significantly affects the inertial sensors.
The velocity estimation with the encoders was most appropriate. But, the positioning result in Figure 12(b) has positioning error because systematic and nonsystematic errors still remained on the less slippery floor.
The result of the proposed system in Figure 12(c) shows reasonable positioning performance. The velocity estimation in Figure 11(e) contains quite large noise. But, the positioning performance of the proposed system is compatible with the positioning with the encoders. Experiments were carried out on a low slip floor, but due to the nature of the odometry, the slip is significantly increased in the acceleration and deceleration section. Therefore, Figure 11(b) seems to be better than Figure 11(e), but it can be seen that the results of Figures 12(b) and 12(c) show almost similar performance when the actual travel distance is measured.
For square motion, the mobile robot moved on a twometerbytwometer square path. The mobile robot moved with onesecond linear acceleration of 1 m/s^{2}, a onesecond constant velocity of 1 m/s, and onesecond linear deceleration at 1 m/s^{2}. The mobile robot turned right 90 degrees at a 90 deg/s angular velocity at each corner. The mobile robot controls velocity by itself with its motor encoders. The experiment for square motion was conducted 10 times. One result with moderate performance was selected as an example. The result of the squaremotion experiment is shown in Figure 13.
The black dotted line is the trajectory of position estimation by the traditional method, which uses an accelerometer, a gyroscope, and a magnetometer. This result is divergent to (in meters) during the experiment’s time of 15 seconds. The purple dashed line is the position estimation from the encoder. This result shows good performance because the experiment was on an even and less slippery floor. The blue line is the trajectory of the position estimation from the proposed method, which uses the ultrasonic anemometer and the IMU. The red point indicates the actual last position of the mobile robot. The position estimation from the encoders and the proposed system was reliable during the squaremotion experiment. The last positions from the encoders and the proposed system are close to the actual last position. Positioning error from the traditional method increases over time. Table 1 shows the stochastic characteristics of each positioning method, aside from the traditional method using IMU because its positioning errors are significantly large. The mean error from odometry is quite large after the wheel of the mobile robot became slippery because of dust. The proposed method shows less mean positioning error and more variance. However, the variance in positioning error is larger than with odometry, and the proposed system measures not wheel revolution but actual movement of the mobile robot. From the experimental results, the ultrasonic anemometer effectively compensates for the accumulated velocity error of the accelerometer. The results show that the proposed system using ultrasonic anemometer and IMU can replace odometry, which has large unpredictable errors in some situations.

5. Conclusion
An IMU offers reasonable data in the short term, but velocity and position errors become large over time due to accumulated errors. In order to reduce error accumulation, odometry or the selfresetting method can be used for dead reckoning. On a slippery or uneven floor, odometry can cause large unpredictable positioning errors. The selfresetting method periodically needs the robot to be stationary to reset its accumulated velocity errors.
In this paper, we proposed a dead reckoning system using an IMU, an ultrasonic anemometer, and a magnetometer. The new approach about the ultrasonic anemometer has been applied. The proposed equation for the ultrasonic anemometer is simpler than the traditional equation. The ultrasonic anemometer prevents velocityerror accumulation in a mobile robot moving through static air. A Kalman filter was designed to combine data from each sensor. The performance of the proposed system was evaluated with experiments. The proposed system showed reasonable performance for velocity and position estimations. The accelerometers reduced noise from the ultrasonic anemometer, and the ultrasonic anemometer prevented accumulation errors in the accelerometers.
Many positioning systems use both absolute and relative positioning systems. The performance of the proposed dead reckoning system can be improved by combining it with a proper absolute positioning method. Then, it could be possible to separate data from the ultrasonic anemometer into the velocity of the mobile robot and the air flow using positioning result.
Conflicts of Interest
The authors declare that there are no conflicts of interest.
Acknowledgments
This work was supported by BK21PLUS, Creative Human Resource Development Program for IT Convergence.
References
 B.S. Cho, W.S. Moon, W.J. Seo, and K.R. Baek, “A dead reckoning localization system for mobile robots using inertial sensors and wheel revolution encoding,” Journal of Mechanical Science and Technology, vol. 25, no. 11, pp. 2907–2917, 2011. View at: Publisher Site  Google Scholar
 J. Borenstein and L. Feng, “Measurement and correction of systematic odometry errors in mobile robots,” IEEE Transactions on Robotics and Automation, vol. 12, no. 6, pp. 869–880, 1996. View at: Publisher Site  Google Scholar
 J. Borenstein and L. Feng, “Gyrodometry: a new method for combining data from gyros and odometry in mobile robots,” in Proceedings of the 13th IEEE International Conference on Robotics and Automation, vol. 1, pp. 423–428, Minneapolis, Minn, USA, April 1996. View at: Google Scholar
 G. A. Piedrahita and D. M. Guayacundo, “Evaluation of accelerometers as inertial navigation system for mobile robots,” in Proceedings of the 3rd IEEE Latin American Robotics Symposium (LARS '06), pp. 84–90, Santiago, Chile, October 2007. View at: Publisher Site  Google Scholar
 M. Ibraheem, “Gyroscopeenhanced dead reckoning localization system for an intelligent walker,” in Proceedings of the International Conference on Information, Networking and Automation (ICINA '10), vol. 1, pp. V167–V172, Kunming, China, October 2010. View at: Publisher Site  Google Scholar
 H.J. Von der Hardt, D. Wolf, and R. Husson, “The dead reckoning localization system of the wheeled mobile robot ROMANE,” in Proceedings of the IEEE/SICE/RSJ International Conference on Multisensor Fusion and Integration for Intelligent Systems, IEEE, December 1996. View at: Publisher Site  Google Scholar
 N. Houshangi and F. Azizi, “Accurate mobile robot position determination using unscented Kalman filter,” in Proceedings of the Canadian Conference on Electrical and Computer Engineering, pp. 846–851, IEEE, May 2005. View at: Publisher Site  Google Scholar
 E. Akeila, Z. Salcic, and A. Swain, “Reducing lowcost INS error accumulation in distance estimation using selfresetting,” IEEE Transactions on Instrumentation and Measurement, vol. 63, no. 1, pp. 177–184, 2014. View at: Publisher Site  Google Scholar
 B. Barshan and H. F. DurrantWhyte, “Inertial navigation systems for mobile robots,” IEEE Transactions on Robotics and Automation, vol. 11, no. 3, pp. 328–342, 1995. View at: Publisher Site  Google Scholar
 J. G. Lee, C. G. Park, and H. W. Park, “Multiposition alignment of strapdown inertial navigation system,” IEEE Transactions on Aerospace and Electronic Systems, vol. 29, no. 4, pp. 1323–1328, 1993. View at: Publisher Site  Google Scholar
 J. Cheng and F. J. De Wan, “A fast initial alignment method for strapdown inertial navigation system on stationary base,” IEEE Transactions on Aerospace and Electronic Systems, vol. 32, no. 4, pp. 1501–1505, 1996. View at: Publisher Site  Google Scholar
 J. Vaganay, M.J. Aldon, and A. Fournier, “Mobile robot attitude estimation by fusion of inertial data,” in Proceedings of the IEEE International Conference on Robotics and Automation, vol. 1, pp. 277–282, IEEE, Atlanta, Ga, USA, May 1993. View at: Publisher Site  Google Scholar
 D. Titterton and J. Weston, Strapdown Inertial Navigation Technology, vol. 17, IET, Stevenage, UK, 2004. View at: Publisher Site
 H. Qi and J. B. Moore, “Direct Kalman filtering approach for GPS/INS integration,” IEEE Transactions on Aerospace and Electronic Systems, vol. 38, no. 2, pp. 687–693, 2002. View at: Publisher Site  Google Scholar
 S. Nassar, Improving the Inertial Navigation System (INS) Error Model for INS and INS/DGPS Applications, Department of Geomatics Engineering, University of Calgary, 2003.
 B.S. Cho, W.J. Seo, W.S. Moon, and K.R. Baek, “Positioning of a mobile robot based on odometry and a new ultrasonic LPS,” International Journal of Control, Automation and Systems, vol. 11, no. 2, pp. 333–345, 2013. View at: Publisher Site  Google Scholar
 N. Kothari, B. Kannan, E. D. Glasgwow, and M. B. Dias, “Robust indoor localization on a commercial smart phone,” Procedia Computer Science, vol. 10, pp. 1114–1120, 2012. View at: Publisher Site  Google Scholar
 A. Widyotriatmo, B. Hong, and K.S. Hong, “Predictive navigation of an autonomous vehicle with nonholonomic and minimum turning radius constraints,” Journal of Mechanical Science and Technology, vol. 23, no. 2, pp. 381–388, 2009. View at: Publisher Site  Google Scholar
 B. Johann and L. Feng, “UMBmark: a benchmark test for measuring odometry errors in mobile robots,” in Mobile Robots X, vol. 2591 of Proceedings of SPIE, International Society for Optics and Photonics, December 1995. View at: Publisher Site  Google Scholar
 H. Surmann, A. Nüchter, and J. Hertzberg, “An autonomous mobile robot with a 3D laser range finder for 3D exploration and digitalization of indoor environments,” Robotics and Autonomous Systems, vol. 45, no. 34, pp. 181–198, 2003. View at: Publisher Site  Google Scholar
Copyright
Copyright © 2017 Woojin Seo and KwangRyul Baek. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.