Flapping insect flight is a complex and beautiful phenomenon that relies on fast, active control mechanisms to counter aerodynamic instability. To directly investigate how freely flying Drosophilamelanogaster control their body pitch angle against such instability, we perturbed them using impulsive mechanical torques and filmed their corrective maneuvers with high-speed video. Combining experimental observations and numerical simulation, we found that flies correct for pitch deflections of up to 40 deg in 29±8 ms by bilaterally modulating their wings' front-most stroke angle in a manner well described by a linear proportional–integral (PI) controller. Flies initiate this corrective process only 10±2 ms after the perturbation onset, indicating that pitch stabilization involves a fast reflex response. Remarkably, flies can also correct for very large-amplitude pitch perturbations – greater than 150 deg – providing a regime in which to probe the limits of the linear-response framework. Together with previous studies regarding yaw and roll control, our results on pitch show that flies' stabilization of each of these body angles is consistent with PI control.

From walking humans to flying insects, many fascinating forms of bio-locomotion are contingent upon robust stabilization control. Implementing this control is particularly difficult in the case of small, flapping-wing insects, as flapping flight is inherently subject to rapidly divergent aerodynamic instabilities (Faruque and Humbert, 2010; Gao et al., 2011; Liu et al., 2010; Pérez-Arancibia et al., 2011; Ristroph et al., 2013; Sun, 2014; Sun et al., 2007; Sun and Xiong, 2005; Taylor and Thomas, 2003; Taylor and Żbikowski, 2005; Xu and Sun, 2013; Zhang and Sun, 2011, 2010). As such, flying insects have evolved stabilization techniques relying on reflexes that are among the fastest in the animal kingdom (Beatus et al., 2015) and robust to the complex environment that insects must navigate (Combes and Dudley, 2009; Dickerson et al., 2012; Ortega-Jimenez et al., 2013; Ravi et al., 2013; Vance et al., 2013).

In particular, pitching instability is a prominent obstacle for flight control in flapping insects. Analytical and numerical modeling suggest that, for many two-winged insects (e.g. flies), periodic flapping couples with longitudinal body motion to produce rapidly growing oscillations of the body pitch angle (Chang and Wang, 2014; Ristroph et al., 2013; Sun, 2014; Sun et al., 2007). This oscillatory instability can be understood as the result of differential drag on the wings due to longitudinal body motion (Ristroph et al., 2013; Sun et al., 2007). For example, if a fly pitches down while hovering, its re-directed lift propels its body forward, causing an increased drag on the wings during the forward stroke relative to the backward stroke. Because the centers of pressure of the fly's wings are always above the body center of mass during normal flapping, this drag asymmetry generates a torque that pitches the fly up. Rather than acting as a restoring torque, the drag – together with the body inertia – pitches the fly up, beyond its initial pitch orientation. The fly then begins to move backwards, and oscillation ensues in the opposite direction. This mechanism results in an undulating instability of the body pitch angle, which doubles over a time scale of ∼9 wingbeats (Sun, 2014). Mitigating the effects of this instability requires flies to actively adjust their wing motion on time scales faster than the growth of these oscillations.

Our work builds upon an already rich corpus of literature on insect flight control, a sizable portion of which addresses the pitch degree of freedom. Experimental studies subjecting tethered insects to both mechanical pitching perturbations and visual pitching stimuli (Dickinson, 1999; Nalbach, 1994; Sherman and Dickinson, 2004, 2003; Taylor and Thomas, 2003; Zanker, 1990) have elucidated stereotyped kinematic responses for pitch correction, including manipulation of wingstroke angle, stroke plane orientation, wingbeat frequency and body configuration. However, tethered insects do not constitute a closed-loop feedback system, as changes to their wing kinematics do not affect their body orientation (Roth et al., 2014). Moreover, in the case of tethered flies, it has been shown that the wing kinematics are qualitatively different from those in free flight (Bender and Dickinson, 2006; Fry et al., 2005). Thus, free-flight studies are necessary for a comprehensive understanding of pitch control. Significant analysis has been performed on freely flying insects executing voluntary maneuvers (Bergou et al., 2010; Ennos, 1989; Fry et al., 2003; Ristroph et al., 2011, 2009) or responding to visual stimuli (Cheng et al., 2011; Muijres et al., 2014; Tammero and Dickinson, 2002; Windsor et al., 2014), but the general challenge of systematically inducing mechanical perturbations on untethered insects has traditionally been a barrier to the study of stabilization reflexes. Some notable exceptions to this include methods of mechanical perturbation using air-flow vortices (Combes and Dudley, 2009; Ortega-Jimenez et al., 2013; Ravi et al., 2013) or gusts of wind (Vance et al., 2013). However, such fluid-impulse methods are difficult to tune, and are thus not ideal for inducing the fast, precise mechanical perturbations that are required for a quantitative understanding of pitch control.

To achieve the necessary speed and precision for measurements of body pitch control, we used a perturbation scheme that has previously been applied to analyzing control of the yaw (Ristroph et al., 2010) and roll (Beatus et al., 2015) degrees of freedom. We mechanically perturbed free-flying Drosophilamelanogaster by gluing small magnetic pins to their dorsal thoracic surfaces and applying short bursts (5–8 ms) of a vertical magnetic field that pitched their bodies up or down. As the flies corrected their orientation, we measured their body and wing kinematics using high-speed video (Fig. 1A).

List of symbols and abbreviations
     
  • CD

    coefficient of drag

  •  
  • CL

    coefficient of lift

  •  
  • CI

    confidence interval

  •  
  • CoM

    center of mass

  •  
  • CoP

    center of pressure

  •  
  • FD

    drag force

  •  
  • FL

    lift force

  •  
  • Ki

    PI controller integral gain

  •  
  • Kp

    PI controller proportional gain

  •  
  • PI

    proportional–integral

  •  
  • second moment of wing area

  •  
  • RMS

    root mean square

  •  
  • S

    wing area

  •  
  • Ut

    wing tip velocity (lab frame)

  •  
  • α

    wing angle of attack

  •  
  • ΔT

    controller latency time

  •  
  • ηw

    wing pitch angle

  •  
  • θb

    body pitch angle

  •  
  • θw

    wing elevation angle

  •  
  • ρ0

    air density

  •  
  • ρb

    body roll angle

  •  
  • τp

    pitch torque

  •  
  • τw

    wing torque

  •  
  • front wingstroke angle

  •  
  • back wingstroke angle

  •  
  • body yaw angle

  •  
  • ϕw

    wingstroke angle

We recorded perturbation events with amplitudes typically ranging from 5 to 40 deg for both pitching up and pitching down. For these perturbations, flies recover 90% of their pitch orientation within ∼30 ms. Moreover, we found that the corrective process is initiated ∼2 wingbeats (∼10 ms) after the onset of the impulsive torque; such a short latency time indicates that this corrective process is a reflexive behavior largely governed by input from the halteres, the flies' rate-gyro-like mechanical sensing organs (Dickinson, 1999; Nalbach, 1994; Pringle, 1948).

To generate corrective pitching torques, flies bilaterally modulate their wings' front-most stroke angle, i.e. they flap their wings more or less in the front to pitch up or down. This corrective mechanism is in general agreement with previous findings on active body pitch stabilization in fruit flies (Chang and Wang, 2014; Dickinson, 1999; Taylor, 2001; Zanker, 1990). We show that flies' modulation of front stroke angle is well modeled by a linear, continuous, proportional–integral (PI) controller with a time delay, ΔT=6±1.7 ms (mean±s.d.). Such simple feedback controller models have been successfully used to describe a wide variety of sensorimotor behavior in many species, including flies (Beatus et al., 2015; Elzinga et al., 2012; Ristroph et al., 2010), hawkmoths (Cheng et al., 2011; Dyhr et al., 2013), electric fish (Cowan and Fortune, 2007; Madhav et al., 2013; Roth et al., 2011), pigeons (Lin et al., 2014) and cockroaches (Cowan et al., 2006; Lee et al., 2008; for comprehensive reviews, see Cowan et al., 2014; Roth et al., 2014; Tytell et al., 2011). Our results indicate that pitch control in fruit flies is an extremely fast and robust process, which can be accurately modeled by a simple controller for a wide range of perturbation amplitudes. Moreover, we found that flies are capable of correcting for pitch deflections of 150 deg or more, a perturbation regime in which the linear controller theory begins to break down. Together with previous results on how flies control yaw (Ristroph et al., 2010) and roll (Beatus et al., 2015), the analysis of pitch control presented here addresses a missing piece in our understanding of how flies control each of their body angles individually.

Fly preparation

We performed each experiment using common fruit flies (D. melanogaster Meigen, females) from an out-bred laboratory stock. Individual flies were anesthetized at 0–4°C, at which point we carefully glued 1.5–2 mm long, 0.15 mm diameter ferromagnetic pins to their notum (dorsal thoracic surface), oriented so that the pin lay in the fly's sagittal plane. The pin is shown in Fig. 1A (false-colored blue) and Fig. 2 (images). Control experiments with untreated flies showed that the addition of the pin did not qualitatively alter flies' flight kinematics. When attached, the pins added ∼20% (∼0.2 mg) to the mass and pitch moment of inertia of the fly, which falls within the range of their natural body mass variations. Moreover, the pin contributed negligibly to off-diagonal components of the flies' inertia tensors, and therefore did not introduce any coupling between the rotational degrees of freedom of the body. The primary effect of the added pin mass was a dorsal shift of ∼0.2 mm in the flies' center of mass, which is accounted for in the calculation of aerodynamic torque. Following Card and Dickinson (2008), the center of mass for untreated flies was assumed to be located at the centroid of the fly body, halfway between the head and tail along the body axis.

Fig. 1.

Pitch perturbation and correction. (A) Snapshots and a 3D model reconstruction from a representative pitch-up event at t=−25.5, 13 and 56.5 ms. The perturbing magnetic field was activated between t=0 and 5.8 ms. For the full video of this event, see Movie 1. The middle snapshot (t=13 ms) corresponds to the maximum pitch-up deflection. The ferromagnetic pin is colored in blue, the fly's long body axis is given by the green arrows, and the red arrow indicates the direction of the perturbation. (B,C) Definitions for the body and wing Euler angle coordinates, respectively. , θb and ρb indicate the body Euler angles (yaw, pitch and roll, respectively), while φw, θw and ηw indicate wing Euler angles (stroke, elevation and pitch, respectively); the stroke plane is indicated by the shaded area (C). Also shown are the lab (B) and body (C) frames of reference. (D–F) Time series of body Euler angles for 18 perturbation events, with the highlighted curves corresponding to the event from A. The yellow bar gives the timing of the magnetic pulse (0–5.8 ms). Body angles are spline-smoothed from raw data, with body yaw in E shifted by its value at t=0.

Fig. 1.

Pitch perturbation and correction. (A) Snapshots and a 3D model reconstruction from a representative pitch-up event at t=−25.5, 13 and 56.5 ms. The perturbing magnetic field was activated between t=0 and 5.8 ms. For the full video of this event, see Movie 1. The middle snapshot (t=13 ms) corresponds to the maximum pitch-up deflection. The ferromagnetic pin is colored in blue, the fly's long body axis is given by the green arrows, and the red arrow indicates the direction of the perturbation. (B,C) Definitions for the body and wing Euler angle coordinates, respectively. , θb and ρb indicate the body Euler angles (yaw, pitch and roll, respectively), while φw, θw and ηw indicate wing Euler angles (stroke, elevation and pitch, respectively); the stroke plane is indicated by the shaded area (C). Also shown are the lab (B) and body (C) frames of reference. (D–F) Time series of body Euler angles for 18 perturbation events, with the highlighted curves corresponding to the event from A. The yellow bar gives the timing of the magnetic pulse (0–5.8 ms). Body angles are spline-smoothed from raw data, with body yaw in E shifted by its value at t=0.

Fig. 2.

Wing kinematics and aerodynamic torques for two representative perturbation events. (A–D) A pitch-up perturbation with maximum amplitude 25 deg; (E–H) a pitch-down perturbation with maximum amplitude 23 deg. Plotted wing kinematics for the left (blue) and right (red) wings as a function of time include: wingstroke angle φw (A,E); wing pitch angle ηw (B,F); and wing elevation angle θw (C,G). Instantaneous aerodynamic torque τ about the flies' pitch axis is also given (D,F). Orange arrows in A and E highlight corrective front strokes, with corresponding arrows in D and H highlighting the changes in pitch torque resulting from the corrective kinematics. Images above A and E show side views of the flies (raw data) at different points during the movie, illustrating the changes in body pitch that accompany a perturbation. White and gray bars indicate forward and back strokes, respectively; the yellow bar corresponds to the perturbation duration (0−5.8 ms).

Fig. 2.

Wing kinematics and aerodynamic torques for two representative perturbation events. (A–D) A pitch-up perturbation with maximum amplitude 25 deg; (E–H) a pitch-down perturbation with maximum amplitude 23 deg. Plotted wing kinematics for the left (blue) and right (red) wings as a function of time include: wingstroke angle φw (A,E); wing pitch angle ηw (B,F); and wing elevation angle θw (C,G). Instantaneous aerodynamic torque τ about the flies' pitch axis is also given (D,F). Orange arrows in A and E highlight corrective front strokes, with corresponding arrows in D and H highlighting the changes in pitch torque resulting from the corrective kinematics. Images above A and E show side views of the flies (raw data) at different points during the movie, illustrating the changes in body pitch that accompany a perturbation. White and gray bars indicate forward and back strokes, respectively; the yellow bar corresponds to the perturbation duration (0−5.8 ms).

Videography and mechanical perturbation

Once 15–30 flies had been prepared as above, we released them into a transparent cubic filming chamber of side length 13 cm. On the top and bottom of the chamber were attached two horizontally oriented Helmholtz coils, which produced a vertical magnetic field. The central region of the chamber was filmed by three orthogonal high-speed cameras (Phantom V7.1) at 8000 frames s−1. The cameras were calibrated using a direct linear transformation scheme detailed in Lourakis and Argyros (2009) and Theriault et al. (2014). When flies entered the filming volume, an optical trigger simultaneously signaled the cameras to record and supplied a 5–8 ms current pulse to the Helmholtz coils. We varied both the strength and duration of the magnetic pulse produced by the coils across experiments. Maximal field strengths reached ∼10−2 T, and most experiments were performed with a pulse that lasted 5.8 ms. The magnetic pulse exerted a torque on the ferromagnetic pin, pitching the fly either up or down.

Of the movies collected from several independent fly batches using the above method, we selected 18 to analyze in full; 16 additional movies were partially analyzed to collect more data on pitch correction time and to observe correction for very large-amplitude perturbations (∼150 deg). We chose movies to fully analyze based on the criteria that (i) the fly was in view of all three cameras for sufficiently long to observe pre- and post-perturbation kinematics, (ii) the perturbation primarily affected the fly's pitch orientation, (iii) the fly was not performing any maneuver other than correction, and (iv) we sampled a wide range of perturbation magnitudes for both pitching up and pitching down across our data set. To glean kinematic data from the raw footage, we used a custom-developed image analysis algorithm detailed elsewhere (Beatus et al., 2015; Ristroph et al., 2009). This 3D hull reconstruction algorithm provided a kinematic description of 12 degrees of freedom for the fly (body orientation and center of mass position, as well as three Euler angles for each wing).

Body and wing kinematics during pitch correction

Representative kinematics for pitch perturbation events are shown in Figs 1 and 2. Fig. 1B shows definitions for the body Euler angle coordinates – yaw (), roll (ρb) and pitch (θb); Fig. 1D–F shows time series of these Euler angles before and after the application of a 5.8 ms magnetic pulse (yellow strip) for 18 perturbation events with at least six different flies. Before the perturbation, flies typically maintained a pitch angle of roughly 50 deg. Perturbations deflected the pitch angle by as much as 40 deg either up or down. While the flies' yaw and roll angles were sometimes altered by the perturbation, pitch was the most consistently and significantly affected degree of freedom immediately following the application of the pulse (at t≈6 ms). Highlighted curves in Fig. 1D–F show an event in which the fly was pitched up by 25 deg, attaining its maximal angular deflection at 15 ms after the onset of the perturbation (Movie 1). By 29 ms it had corrected for 90% of the pitch deflection. The maximum pitch velocity from the perturbation was 2400 deg s−1.

The wing kinematics for two representative perturbations, one in which the fly was pitched up by 25 deg (highlighted in Fig. 1D–F) and another in which the fly was pitched down by 23 deg, are shown in Fig. 2. Wing Euler angles – stroke (φw), pitch (ηw) and elevation (θw) – are defined in Fig. 1C. In general, the wing kinematics we observed following the perturbation were left/right symmetric. For the pitch-up event, ∼10 ms, or 2 wingbeats, after the onset of the perturbation, the minima of the wingstroke angles shifted upward for both the left and right wing (Fig. 2A, orange arrows). During a given wingbeat cycle, the minimum of the stroke angle for each wing corresponds to its front-most position. We refer to the average of the front-most positions for the left and right wings as the front wingstroke angle, . By t=15 ms, the fly in Fig. 2A–D had increased its from its pre-perturbation value by ∼25 deg. Physically, this means that the fly was significantly reducing the amplitude of its ventral stroke, i.e. flapping less forward; the duration of this increase in was 3 wingbeats. We did not observe any shifts in the fly's back-most stroke angle, , during the correction maneuvers.

Conversely, for the pitch-down event in Fig. 2E–H, ∼10 ms after the perturbation onset, the pitched-down fly began to decrease its . This corresponds to the fly increasing the amplitude of its ventral stroke, i.e. flapping further forward. Again, there appeared to be little to no change in the back stroke angle.

Aerodynamic forces and torques

Intuitively, the relationship between front stroke angle and pitching torque can be understood as follows. To within a good approximation, the net aerodynamic force generated by a flapping wing is directed perpendicular to the wing's surface (Dickinson et al., 1999; Sane and Dickinson, 2001), so that portions of the wingstroke during which the wing is in the front half of the stroke plane (φw≤90 deg) generate pitch-up torques, while portions in the back half (φw≥90 deg) generate pitch-down torques (Fig. 2D,H, Fig. 3). During non-maneuvering flight, these torques cancel over a wingstroke. By biasing a wingstroke so that it spends a smaller fraction of the stroke period in the forward position, a smaller pitch-up torque is generated during that cycle, such that the net pitch torque will be directed downward. Conversely, by increasing the front stroke angle, and thus increasing the portion of the stroke spent in the front position, flies can generate a net pitch-up torque over the course of a full stroke. This can be observed in Fig. 2 (orange arrows) and Fig. 4, in which active adjustments to front stroke angle resulted in net corrective pitching torques over the course of individual wingbeats.

Fig. 3.

Schematic diagram of moment arm and aerodynamic force direction for wing-generated torques. Left, a fly with its wings in the back half of the stroke plane generating a pitch-down torque (red arrow) and, right, the wings in the front half of the stroke plane, generating a pitch-up torque (red arrow). Moment arms, r (green arrows), are defined as the vector from the fly's body center of mass (CoM) to the wing's center of pressure (CoP); our data showed that the CoP is consistently above the CoM during flight. The net aerodynamic force vector, Fw (blue arrows), is oriented normal to the wing surface; however, the relevant components of Fw for pitching torques lie in the sagittal plane (as in the diagram). Moreover, for symmetric flapping, the components of Fw out of the sagittal plane cancel between the left and right wings.

Fig. 3.

Schematic diagram of moment arm and aerodynamic force direction for wing-generated torques. Left, a fly with its wings in the back half of the stroke plane generating a pitch-down torque (red arrow) and, right, the wings in the front half of the stroke plane, generating a pitch-up torque (red arrow). Moment arms, r (green arrows), are defined as the vector from the fly's body center of mass (CoM) to the wing's center of pressure (CoP); our data showed that the CoP is consistently above the CoM during flight. The net aerodynamic force vector, Fw (blue arrows), is oriented normal to the wing surface; however, the relevant components of Fw for pitching torques lie in the sagittal plane (as in the diagram). Moreover, for symmetric flapping, the components of Fw out of the sagittal plane cancel between the left and right wings.

Fig. 4.

Mechanism for generating corrective pitch torques. (A,B) Data correspond to the pitch-up event shown in Fig. 1D–F and Fig. 2A–D; (C,D) data correspond to the pitch-down event shown in Fig. 2E–H. (A,C) Wingstroke-averaged pitch torque (black line, spline smoothed) and negative front stroke angle deviation (orange dots) versus time. Front stroke angle deviation is defined by . We flip the sign of to illustrate the correlation with time-averaged body torque (τp). (B,D) Body pitch angle (Δθb, red line) versus time. As in previous plots, the yellow strip in A–D corresponds to the duration of the magnetic pulse. Images above A and C are raw data snapshots from the overhead camera of the flies facing to the right, showing the front-most stroke angle for the indicated time points.

Fig. 4.

Mechanism for generating corrective pitch torques. (A,B) Data correspond to the pitch-up event shown in Fig. 1D–F and Fig. 2A–D; (C,D) data correspond to the pitch-down event shown in Fig. 2E–H. (A,C) Wingstroke-averaged pitch torque (black line, spline smoothed) and negative front stroke angle deviation (orange dots) versus time. Front stroke angle deviation is defined by . We flip the sign of to illustrate the correlation with time-averaged body torque (τp). (B,D) Body pitch angle (Δθb, red line) versus time. As in previous plots, the yellow strip in A–D corresponds to the duration of the magnetic pulse. Images above A and C are raw data snapshots from the overhead camera of the flies facing to the right, showing the front-most stroke angle for the indicated time points.

To quantify the effect of changing the front stroke angle, we calculated the pitching torque generated by the wings during the maneuvers shown in Fig. 2 using the full 3D fly kinematics. To calculate the aerodynamic force generated by the wings, we used a quasi-steady aerodynamic force model that was previously calibrated on a mechanical, scaled-up fly model (Dickinson et al., 1999; Sane and Dickinson, 2001). This model gives the lift (FL) and drag (FD) forces generated by the wings as:
(1)
(2)

where CL and CD are the wing's lift and drag coefficients, respectively, and are given as functions of wing angle of attack (α) by Sane and Dickinson (2001); S is the wing area; is the non-dimensional radius of the second moment of wing area (given as 0.313 by Cheng et al., 2009); ρ0 is the density of air; and Ut is the wingtip velocity as measured in the lab frame. Drag is directed anti-parallel to the wing tip velocity, and lift is perpendicular to both drag and the wing span vector. The total aerodynamic force is the vector sum of the lift and drag forces. While this is a simple method for calculating forces on flapping wings, we found that it quantitatively captures the relevant force production for both pre- and post-perturbation wing kinematics. We tested the effect of adding rotational forces to the aerodynamic model (Sane and Dickinson, 2001), which should give the next largest contribution to force and account for the force peaks at wingstroke reversal, but found negligible changes to our calculation results.

From these lift and drag forces, we calculated aerodynamic torques exerted by the wings on the body (Fig. 2D,H). The moment arm for the torque is given by the vector from the fly's center of mass to the wing's center of pressure, assumed to be in the chord center, 70% along the length of the wing's span (Fig. 3) (Cheng et al., 2009; Fry et al., 2005). During the active correction of the pitch-up perturbation, the fly's wings generated less upward torque, resulting in a net downward bias for the torque on the body (orange arrows in Fig. 2D). Similarly, the pitched-down fly in Fig. 2H generated net upward pitching torque during active correction (orange arrows in Fig. 2H). The corrective torques for these two events are well correlated with the measured modulations of front stroke angle (Fig. 4A,C), a trend that we observed across all perturbation events (see Results, ‘The corrective effect of stroke angle over a range of perturbations’).

Interestingly, after the flies generated a corrective torque for 2–3 wingbeats, we also observed a few wingbeats in which they generated net torque in the opposite direction (Fig. 2D,H, Fig. 4A,C). As with the initial corrective torque, this subsequent counter-torque arose from modulations of the front stroke angle, evident in Fig. 4A,C. The counter-torque acted to brake the corrective pitching motion, mitigating the overshoot in body pitch angle caused by the initial correction response. This allows for faster correction times, as the initial corrective maneuver can generate larger torques, and thus more quickly return the fly to pitch angles near its original orientation. In movies that allowed us to track the fly for long periods after the perturbation, we observed that the front stroke angle and the net aerodynamic torque often oscillated with decaying amplitude and a period of ∼3–4 wingbeats.

Importantly, passive damping of pitch motion contributed negligibly to the correction maneuvers we observed, in contrast to the case of yaw stability (Cheng et al., 2009; Hedrick et al., 2009; Hesselberg and Lehmann, 2007; Sun, 2014; Warrick et al., 2012). The characteristic time scale at which passive pitch damping becomes significant has been estimated as ∼80 ms (Ristroph et al., 2013), much longer than the entire correction maneuver. Taken together, our results indicate that pitch correction for flies is an active process involving modulation of .

The importance of stroke angle relative to other degrees of freedom

To assess the effect of front stroke angle modulation on body pitch correction, we calculated the changes to body pitch angle, Δθb, generated by changes in wing kinematics. We isolated the corrective effect of each of the three wing kinematic variables by first identifying wingstrokes that correspond to both non-maneuvering (no net torque) and corrective flight. We calculated the Δθb resulting from all eight combinations of corrective (red) and non-corrective (blue) kinematics (color code in Fig. 5). For example, a triplet of squares with color combination red–blue–blue corresponds to a wingbeat with the stroke angle taken from a maneuvering wingbeat, and both wing pitch and elevation angles taken from a non-maneuvering wingbeat; its coordinate along the horizontal axis gives the calculated change in body pitch angle resulting from such a wingbeat. Our calculation assumed rigid wings attached to a stationary body with the geometry as in Fig. 3, and we determined the net change in body pitch angle using numerical integration and assuming an initial condition with zero pitch angular velocity.

Fig. 5.

The relative importance of wing kinematic variables. The change in body pitch angle, Δθb, resulting from different combinations of wing Euler angles for both a pitch-up and a pitch-down event. Each triplet of squares corresponds to a unique combination of wing angle kinematics taken from either maneuvering (red) or non-maneuvering (blue) wingbeats. These combinations of Euler angle kinematics are used with our quasi-steady aerodynamic model to calculate pitching torques and pitch angle deflection over a single wingbeat. Wing kinematics φw (stroke), ηw (wing pitch) and θw (elevation) correspond to individual squares within each rectangular triplet. The wing angles are defined in Fig. 1C. The top plot corresponds to kinematics taken from the pitch-up event in Figs 2 and 4; the bottom plot corresponds to kinematics taken from the pitch-down event in Figs 2 and 4. A color key and a case example for a particular combination of stroke parameters (red–blue–blue) are shown at the bottom.

Fig. 5.

The relative importance of wing kinematic variables. The change in body pitch angle, Δθb, resulting from different combinations of wing Euler angles for both a pitch-up and a pitch-down event. Each triplet of squares corresponds to a unique combination of wing angle kinematics taken from either maneuvering (red) or non-maneuvering (blue) wingbeats. These combinations of Euler angle kinematics are used with our quasi-steady aerodynamic model to calculate pitching torques and pitch angle deflection over a single wingbeat. Wing kinematics φw (stroke), ηw (wing pitch) and θw (elevation) correspond to individual squares within each rectangular triplet. The wing angles are defined in Fig. 1C. The top plot corresponds to kinematics taken from the pitch-up event in Figs 2 and 4; the bottom plot corresponds to kinematics taken from the pitch-down event in Figs 2 and 4. A color key and a case example for a particular combination of stroke parameters (red–blue–blue) are shown at the bottom.

We performed this analysis for data from two different perturbation events: the pitch-up and pitch-down events in Figs 2 and 4. The grouping of the points in Fig. 5 indicates that body pitch correction is most closely associated with changes to wingstroke angle. The red–blue–blue point, corresponding to corrective stroke angle but non-maneuvering wing pitch and elevation, achieves at least 60% of the body pitch correction, consistent with Muijres et al. (2014). Moreover, the only combinations that account for more than 30% of the total correction include the stroke angle from a maneuvering wingbeat (points of the form red–xx). These results motivate a minimal model for body pitch stabilization that considers only variations in front stroke angle to drive pitch correction.

The corrective effect of stroke angle over a range of perturbations

To further flesh out the relationship between front stroke angle and corrective torque, we plotted the maximum measured corrective pitch acceleration generated by the fly in each of the 18 maneuvers as a function of the corresponding change in front stroke angle measured at that time, (Fig. 6A). The maximum pitch acceleration was measured at the extremum of θb, using a quadratic polynomial fit. The plot demonstrates a strong correlation between changes in front stroke angle and corrective acceleration (linear R2=0.87). Across our data set of 18 perturbation events, flies increased (flap less forward) to pitch themselves down, and decreased (flap further forward) to pitch themselves up.

Fig. 6.

Statistics from 18 perturbation events. (A,B) Maximum corrective pitch acceleration () generated by the fly as a function of change in front stroke angle () and back stroke angle (). The gray line in A is the calculated pitch acceleration for simplified wing kinematics. This calculation has no fit parameters, and is based on morphological parameters (Cheng et al., 2009), wing kinematics (Chang and Wang, 2014) and the quasi-steady aerodynamic model. (C) A histogram of latency times across our data set, with latency time defined as the time between the onset of the magnetic perturbation and the beginning of a measurable corrective wing response (±4 deg change in front stroke angle). The orange background corresponds to the mean delay time (±s.d.) obtained from our controller model fits. (D) The time for pitch correction – defined as the time it takes for the fly to recover 90% of its original pitch orientation – plotted as a function of the maximum pitch deflection (Δθb,max) for each perturbation event. Solid and dashed lines give the mean±s.d. The lack of discernible correlation (R2=0.0093) is consistent with linear control.

Fig. 6.

Statistics from 18 perturbation events. (A,B) Maximum corrective pitch acceleration () generated by the fly as a function of change in front stroke angle () and back stroke angle (). The gray line in A is the calculated pitch acceleration for simplified wing kinematics. This calculation has no fit parameters, and is based on morphological parameters (Cheng et al., 2009), wing kinematics (Chang and Wang, 2014) and the quasi-steady aerodynamic model. (C) A histogram of latency times across our data set, with latency time defined as the time between the onset of the magnetic perturbation and the beginning of a measurable corrective wing response (±4 deg change in front stroke angle). The orange background corresponds to the mean delay time (±s.d.) obtained from our controller model fits. (D) The time for pitch correction – defined as the time it takes for the fly to recover 90% of its original pitch orientation – plotted as a function of the maximum pitch deflection (Δθb,max) for each perturbation event. Solid and dashed lines give the mean±s.d. The lack of discernible correlation (R2=0.0093) is consistent with linear control.

The correlation between and corrective pitch acceleration in Fig. 6A is also predicted by a calculation based on the quasi-steady aerodynamic force model in Eqns 1 and 2 (Fig. 6A, gray line). To calculate the aerodynamic pitch torques, we used a simplified wing kinematic model similar to that in Chang and Wang (2014) in which only the front stroke angle is varied (see Appendix). We averaged the computed torques over a wingbeat, and divided by the moment of inertia to obtain pitch acceleration. The calculation relied only on the wing kinematics and fly morphology (Cheng et al., 2009), and had no fitted parameters. The calculation quantitatively reproduced the measured pitch acceleration, further corroborating a model for pitch control that includes only modulation of .

To rule out an alternative corrective mechanism, based on modulation of back stroke angle, we plotted corrective pitch acceleration as a function of (Fig. 6B) and found no discernible correlation between these two variables.

Correction time scales

We also analyzed the pitch correction time scales across our data set in terms of both the response latency and the overall duration of the correction maneuver. Fig. 6C shows a histogram of latency times for each corrective maneuver, defined as the time between the onset of the perturbation and the first measurable change in the front stroke angle ( deg) for all 18 perturbation events. The mean latency time was 9.9±2.1 ms (mean±s.d., N=18), corresponding to ∼2 wingbeats. Fig. 6D plots the total correction time for each maneuver as a function of the maximum body pitch deflection in each perturbation event. We define the correction time as the time between the onset of the perturbation and the fly correcting 90% of the pitch deflection. The mean correction time was 29±8 ms (mean±s.d., N=32). Finally, we found that the correction time is weakly correlated with the perturbation amplitude (linear R2=0.093), which is consistent with a linear control model.

Control-theory model

We used a control-theoretic framework to describe the flies' strategy for pitch stabilization. In particular, we modeled actuated changes to the front stroke angle as the output of a PI controller with time delay ΔT, for which the input is body pitch velocity (block diagram in Fig. 7A). The response is given by:
(3)
Fig. 7.

Control theory model. (A) A block diagram for the simplified proportional–integral (PI) controller model. The ‘Fly’ box denotes the fly dynamics, which were determined by the Newton–Euler equations (a simplified model for these dynamics is detailed in the Appendix). The halteres sense the body's pitch angular velocity, . This pitch velocity signal is subject to a time delay, ΔT, and split into two branches. One branch is multiplied by Ki and integrated to yield pitch displacement, while the other branch is multiplied by Kp. These signals are recombined as an output, , that adjusts the front stroke angle of the wings and results in a corrective wing torque (τw). The corrective wing torque affects the fly's motion and closes the feedback loop in the stabilization controller. (B) Measured front stroke angle as a function of time for the pitch-up event in Figs 1, 2 and 4 (orange dots) compared with the output of the fitted PI controller model (blue line). The relative contributions from the proportional (P) and integral (I) terms are shown by the gray solid line and the brown dashed line, respectively. Confidence intervals (shaded blue region) were calculated based on a χ2 test for the fitting residuals in the control parameter space.

Fig. 7.

Control theory model. (A) A block diagram for the simplified proportional–integral (PI) controller model. The ‘Fly’ box denotes the fly dynamics, which were determined by the Newton–Euler equations (a simplified model for these dynamics is detailed in the Appendix). The halteres sense the body's pitch angular velocity, . This pitch velocity signal is subject to a time delay, ΔT, and split into two branches. One branch is multiplied by Ki and integrated to yield pitch displacement, while the other branch is multiplied by Kp. These signals are recombined as an output, , that adjusts the front stroke angle of the wings and results in a corrective wing torque (τw). The corrective wing torque affects the fly's motion and closes the feedback loop in the stabilization controller. (B) Measured front stroke angle as a function of time for the pitch-up event in Figs 1, 2 and 4 (orange dots) compared with the output of the fitted PI controller model (blue line). The relative contributions from the proportional (P) and integral (I) terms are shown by the gray solid line and the brown dashed line, respectively. Confidence intervals (shaded blue region) were calculated based on a χ2 test for the fitting residuals in the control parameter space.

Eqn 3 states that adjustment of the front stroke angle () at a given time t is given by a linear combination of the body pitch angle deviation from the pre-perturbation orientation (Δθb) and body pitch velocity (), both measured at an earlier time t−ΔT. The parameters Kp and Ki are the proportionality constants that determine the relative weights of body pitch angle and pitch velocity. Note that the same controller could be termed a proportional–derivative (PD) controller if the input to the controller was chosen to be the body pitch angle. Because the fly halteres are known to measure body angular velocities (Dickinson, 1999; Nalbach, 1994; Pringle, 1948), we choose the PI nomenclature. We did not consider controller models that depend on angular acceleration (proportional–integral–derivative, PID) based on previous studies that have shown flies' corrective pitch response to be insensitive to angular acceleration (Dickinson, 1999).

Using measured values for, Δθb and , we fitted for the parameters Kp, Ki and ΔT. The three parameters were fitted for each movie individually, using one data point per wingstroke. The fit was performed by evaluating Eqn 3 on a dense 3D grid in parameter space and finding the global minimum for the sum of squared residuals between the control model and real data. The results of a controller fit for the pitch-up event in Figs 1, 2 and 4 are shown in Fig. 7B (orange dots, blue line). We found excellent quantitative agreement between the controller fit and our measured data with a root mean square (RMS) error of 1.9 deg, which is of the order of the measurement uncertainty for (Ristroph et al., 2009). The controller model captures not only the sharp rise in in response to the perturbation but also the subsequent decrease in corresponding to the braking counter-torque that slows the fly's downward pitching motion (see Results, ‘Aerodynamic forces and torques’). The fast rise time of the response can be attributed to the proportional term.

We applied the same fitting process to nine movies. We found the values of fitted control parameters (Table 1) to be Ki=0.3±0.15, Kp=7±2.1 ms and ΔT=6±1.7 ms (means±s.d.), with an average RMS fitting error of 3.0 deg. The mean value of ΔT corresponds to ∼1 wingbeat. Fig. 6C shows the region corresponding to mean ΔT±1 s.d. (highlighted in orange) compared with measured latency times. The confidence interval (CI) size is large relative to the fitted control parameters (>50% in some cases). The large confidence intervals, combined with the accuracy of the fit, indicate that the controller output is insensitive to the choice of model parameters, i.e. that the controller gains and time delay do not require fine tuning.

Table 1.

Fit results for proportional–integral (PI) controller model with confidence intervals (CI) for each parameter

Fit results for proportional–integral (PI) controller model with confidence intervals (CI) for each parameter
Fit results for proportional–integral (PI) controller model with confidence intervals (CI) for each parameter

Numerical simulation

To corroborate our experimental evidence for the PI controller, we performed a dynamical simulation of a mechanically perturbed fruit fly. The simulation solved the equations of motion for the pitch, longitudinal and vertical degrees of freedom, assuming the fly's geometry, simplified wing kinematics and the quasi-steady aerodynamic force model detailed above (for details, see Appendix). The body pitch angle over time for simulated flies implementing different control strategies is shown in Fig. 8. The four simulated control schemes shown are: (i) proportional–integral (PI, blue), (ii) proportional (P, green), (iii) integral (I, orange) and (iv) no control (red). To determine parameters for the simulated controllers, we fitted the parameters of each model to experimental data. To mimic experimental conditions, we imposed a 5 ms external mechanical torque on the simulated flies (yellow strip), with magnitude comparable to our real system.

Fig. 8.

Numerical simulation results for different controller models. (A) Time series of body pitch angle θb for simulated flies implementing proportional–integral (PI, blue) and proportional (P, green) control. (B) Time series of body pitch angle for simulated flies implementing integral (I, orange) and no (None, red) control. The time axes for A and B are the same, but the range of the pitch angle axis differs significantly between the two. The gray region in A and B corresponds to 45±2 deg, where 45 deg is the reference body pitch angle for each controller. Hence, curves returning to and remaining within this region indicate successful control. The yellow strip indicates the duration of the mechanical perturbation (0–5 ms).

Fig. 8.

Numerical simulation results for different controller models. (A) Time series of body pitch angle θb for simulated flies implementing proportional–integral (PI, blue) and proportional (P, green) control. (B) Time series of body pitch angle for simulated flies implementing integral (I, orange) and no (None, red) control. The time axes for A and B are the same, but the range of the pitch angle axis differs significantly between the two. The gray region in A and B corresponds to 45±2 deg, where 45 deg is the reference body pitch angle for each controller. Hence, curves returning to and remaining within this region indicate successful control. The yellow strip indicates the duration of the mechanical perturbation (0–5 ms).

Among the four tested models, PI control is the only one that is consistent with the fast, robust pitch control that we observed experimentally (Fig. 8A). We found that flies with no control or I control exhibited large, rapid oscillations of body pitch angle (Fig. 8B), while flies with P control exhibited slightly smaller, long time scale oscillations (Fig. 8A). With I, P and no control, the simulated fly failed to remain aloft and rapidly lost altitude. In contrast, simulated flies implementing PI control corrected their orientation over time scales similar to those in our experimental data, maintained pitch stability over long times, and remained aloft. The general features of each control scheme showed little sensitivity to the values of the control parameters, in agreement with the large confidence intervals of the fitted control parameters obtained above (see ‘Control-theory model’).

Extreme perturbations

In addition to the 18 perturbation events analyzed in full (Figs 1, 6), we examined two large-amplitude perturbation events of two additional flies. Snapshots and time courses of body pitch angle are shown in Fig. 9 for a pitch-up and a pitch-down event, both with maximum pitch deflection greater than 130 deg (Movie 2). Remarkably, both flies performed successful correction maneuvers, although they were not in-frame long enough to observe them returning to their original orientation. The correction time for both large-amplitude events (>50 ms) was longer than the correction times shown in Fig. 6, which can be attributed to the fact that the controlled quantity is biologically constrained: front stroke angle is limited, for instance, by the angle at which the body or the other wing obstructs a wing's forward motion. When we input the body pitch kinematics for the events in Fig. 9 into our PI controller model, the model predicted changes to in excess of 100 deg – a value that is physiologically impossible in the forward direction and not observed in the backward direction. Assuming the flies' corrective response is bounded by  deg, a physiologically reasonable estimate, our numerical simulation predicted a response time of ∼70 ms for a perturbation Δθb=−150 deg, in excellent agreement with the experimental data.

Fig. 9.

Large perturbations. Overlaid snapshots from raw data and time series of body pitch angle for (A) pitch-up perturbation and (B) pitch-down perturbation (Movie 2). The pitched-up fly reaches a maximum pitch deflection of 130 deg at ∼20 ms after the onset of the perturbation. The pitched-down fly reaches a maximum pitch deflection of −155 deg after ∼30 ms. The loss of altitude during the correction is evident in both cases and shown to scale. Yellow strips indicate the 5.8 ms magnetic pulse; red arrows indicate the direction of each perturbation.

Fig. 9.

Large perturbations. Overlaid snapshots from raw data and time series of body pitch angle for (A) pitch-up perturbation and (B) pitch-down perturbation (Movie 2). The pitched-up fly reaches a maximum pitch deflection of 130 deg at ∼20 ms after the onset of the perturbation. The pitched-down fly reaches a maximum pitch deflection of −155 deg after ∼30 ms. The loss of altitude during the correction is evident in both cases and shown to scale. Yellow strips indicate the 5.8 ms magnetic pulse; red arrows indicate the direction of each perturbation.

Front stroke angle as the controlled quantity

We showed that front stroke angle modulation is the primary mechanism for body pitch control in fruit flies, consistent with previous experiments (Dickinson, 1999; Ristroph et al., 2013; Taylor, 2001; Zanker, 1990), and in the same spirit as other proposed mechanisms that include modulation of the mid-stroke angle (Chang and Wang, 2014). Calculated aerodynamic forces predicted that changes to the back stroke angle could produce corrective pitching torques in the same way that changes to front stroke angle do; the fact that we did not observe this in the data hints that morphological constraints favor modulation of the front stroke angle. Our computational results (Figs 6, 8) showed that a minimal model, which only incorporates changes to front stroke angle, and uses control parameters extracted from fits to our measurements, is the simplest linear, continuous model capable of stabilizing the body pitch angle on time scales similar to those observed in the experiments.

Kinematic variables other than may also contribute to pitch correction. Previous studies have associated changes to stroke plane elevation (Zanker, 1988), wingbeat frequency (Dickinson, 1999) and body posture (Taylor, 2001) with pitch correction. In particular, we observed transient alterations in both wing pitch and elevation angle during corrective maneuvers. Fig. 5 suggests that, when combined with modulation of front stroke angle, changes to wing pitch and elevation angle can account for up to 40% of body pitch correction, consistent with Muijres et al. (2014). The detailed role of these kinematic variables in pitch control and whether they are actively or passively actuated remains unknown (Bergou et al., 2010).

Long-term stabilization

Importantly, the PI controller model presented here, which assumes a ‘dead reckoning’ method of angle estimation based on angular velocity input, cannot account for pitch control on long time scales because of integration and sensor errors affecting the estimation of the absolute pitch angle. Long-term pitch control in this framework requires direct measurement of the pitch angle, as could be achieved by the visual system at longer time scales (see Dyhr et al., 2013 for a discussion of time delays in the visual response system). An interplay between the haltere and visual systems (as in Huston and Krapp, 2009; Sherman and Dickinson, 2004) would thus be necessary for comprehensive pitch stability. An intriguing alternative to the dead reckoning assumption is that the fly could implement a model-based estimator to make measurements of its absolute pitch angle using only information about its pitch velocity (N. Cowan, personal communication). Such a model-based estimator is indeed possible: the linearized, flying insect pitch dynamics in Ristroph et al. (2013) give rise to an observability matrix that is full rank, indicating that body pitch angle should be observable with only angular velocity feedback. Whether this method for control is actually used by flies would be the subject of further research. The PI model presented here, however, can completely and accurately account for the fly's fast reflex response, which stabilizes it against rapid pitch perturbations.

Discrete versus continuous control models

The periodic motion of wing flapping introduces inherent discreteness to insect flight. For processes occurring on time scales comparable to a wingstroke period – like the perturbations and maneuvers we reported here – we expect discrete effects to be more pronounced. In particular, modulations of front stroke angle can, by definition, occur only once per wingbeat. Because perturbations can be induced at any time during the wingbeat, but the actuated kinematics are discretely constrained, latency times for correction depend on the phase of the perturbation relative to the wingstroke. Latency times will be bounded from below by the flies' neural response time, but could potentially be as much as 1 wingbeat longer as a result of the phase of the perturbation within the wingbeat.

Measured latency times can also depend on discrete sensing. The temporal sampling resolution with which flies can measure mechanical perturbations is likely determined by motion of their halteres, the rate-gyro sensory organs used in the fast perturbation response (Dickinson, 1999). Dipteran halteres beat at the wing frequency and use Coriolis forces to measure body angular velocities (Nalbach, 1994; Pringle, 1949). The largest sensitivity to mechanical perturbations is likely to occur at times during the fly's mid-stroke, when Coriolis forces on the halteres are the largest (Nalbach, 1994). Sensing at discrete times introduces a second relevant phase for correction latency time: the phase of the perturbation relative to sensing. Similar to discrete actuation, discrete sensing would lead to latency times longer than the neural response time. Moreover, even during the fly's mid-stroke, its halteres only have finite sensitivity. It is likely that there exists some threshold for angular velocities that are large enough to elicit a control response (Fox and Daniel, 2008).

The continuous PI controller model does not account for the aforementioned effects of discrete actuation, discrete sampling or sensing threshold. Hence, the measured latency time should constitute an upper bound for the delay time that we obtain from the controller model. Indeed, the time delay from our controller model (ΔT=6±1.7 ms) is roughly 1 wingbeat shorter than the measured latency time (9.9±2.1 ms).

Despite the inherent discreteness of the fly control systems, the continuous PI controller model quantitatively captures the behavior of flies in response to pitch perturbations (Fig. 7). This quantitative agreement leads to an interesting open question: under what conditions does it become necessary to use a discrete controller model to describe flight stabilization? Previous studies on the legged locomotion of cockroaches have shown quantitative consistency between discrete and continuous control models in the context of wall following (Lee et al., 2008). To address this question in the context of flight stabilization would require precise perturbation timing, in order to probe the short time scales at which discretization becomes relevant. Such an analysis could provide significant insight into the timing and thresholding of fruit fly reflexes.

Physiological basis for pitch control

Both the mechanism and timing of the pitch correction indicate a likely candidate muscle for control actuation: the first basalare muscle (b1), as suggested by previous studies (Chang and Wang, 2014; Fayyazuddin and Dickinson, 1999). Among dipteran flight control muscles, b1 is unique in that it is active during every wingstroke (Heide, 1983; Heide and Götz, 1996; Miyan and Ewing, 1985), which would allow for the wingbeat time scale pitch control that we observed. Moreover, b1 activity is strongly correlated with modulations of ventral stroke amplitude, i.e. changes in (Dickinson and Tu, 1997; Walker et al., 2014). In blowflies, b1 activity is also correlated with changes in elevation angle during the ventral stroke (Balint and Dickinson, 2004), which could explain the slight shifts in elevation angle that we observed during correction (Fig. 2B,F). Our results indirectly support the hypothesis that the b1 muscle is responsible for pitch control through the regulation of . Testing flies with disabled or altered b1 muscles could provide an avenue for confirming the role of b1 in the body pitch control process.

Linear control of body orientation

In addition to the results on fruit flies reported here, PI control has also been identified in pitch control for hawkmoths (Cheng et al., 2011; Windsor et al., 2014). The anatomical similarities found across species suggest that pitch instability is an obstacle faced by many flapping insects (Ristroph et al., 2013; Sun, 2014; Sun et al., 2007); a natural question raised by these collective findings is whether PI control is a generic feature of pitch stabilization in flying insects. Beyond flying insects, what we refer to as PI control has also been observed in fast obstacle avoidance in pigeons (Lin et al., 2014) and antenna-based wall following in cockroaches (Cowan et al., 2006; Lee et al., 2008). Future research on the ubiquity of PI control could have fascinating implications for the evolution of flight stabilization and sensorimotor control mechanisms.

In the larger context of comprehensive flight stabilization, our results, together with previous studies on yaw (Ristroph et al., 2010) and roll (Beatus et al., 2015) control in fruit flies, show that the strategies flies use to control each of their body Euler angles can be modeled as PI controllers. However, the overarching structure in which these three individual controllers are embedded is still unknown. Given the non-commutativity of rotations in 3D, the relationship between controllers that measure different angular coordinates is likely to be non-trivial. For example, a previous study has shown two cases where the fly's response to certain perturbation could not be explained by a superposition of linear controllers (Beatus et al., 2015). First, in response to perturbations that simultaneously affected both the roll and yaw angles, preferential correction for roll over yaw was observed, hinting that the control of these two degrees of freedom is coupled. Second, it was shown that the response of flies to extreme perturbations consisting of multiple rotations along roll cannot be explained by a linear control model. Taken together, these results hint at a complex and intriguing control architecture used by flies to stabilize their orientation that may depend on the amplitude and timing of the perturbations along each axis.

Finally, an understanding of the relationship between control of different Euler angles could have profound implications for how the fly encodes information about its body orientation. In the case of vision, organism-specific demands have spurred the development of novel, specialized neural structures in both mammals (Hafting et al., 2005; Yartsev et al., 2011) and insects (Ofstad et al., 2011; Seelig and Jayaraman, 2013). Pioneering work on information processing from halteres has suggested similar morphology/function relationships for the gyroscopic rate sensing in insects (Fox et al., 2010). Connecting such analyses with the resultant control structure observed in free-flight behavioral experiments could provide a window into the most basic ways in which flies sense and interpret the world.

Appendix

Simplified wing kinematics

For both the calculation of aerodynamic torques in Fig. 6A and the dynamical simulation in Fig. 8, we use an analytic form for simplified wing kinematics taken from Chang and Wang (2014). These kinematics closely resemble the motion of real fly wings, but are simple enough to write down concisely:
(A 1)
(A 2)
(A 3)

The wing Euler angles are defined in Fig. 1C. The terms in Eqns A1–3 are defined as follows: φ0, η0 and θ0 are angle offsets; φm, ηm and θm are amplitudes; K and C are waveform parameters – K tunes the stroke angle from pure sine wave to triangle wave, while C tunes the wing pitch angle from sinusoid to square wave; ω is the wingbeat frequency; and δη and δθ are phase offsets.

For the parameter values used in the calculation in Fig. 6A, see Table 2. Note that in the main text we refer to front and back stroke angles ( and ), which are related to φm and φ0 by the linear relationships:
(A 4)
(A 5)
Table 2.

Model parameter symbols, definitions and values

Model parameter symbols, definitions and values
Model parameter symbols, definitions and values

Numerical simulation

Using the simplified wing kinematics above and the quasi-steady aerodynamic model detailed in the main text (Eqns 1 and 2), our numerical simulation solves the Newton–Euler equations for vertical (z) and forward (x) motion of the body center of mass, as well as pitch rotational (θb) motion:
(A 6)
(A 7)
(A 8)

where M is the fly body mass, Ipitch is its pitch moment of inertia, g is the gravitational acceleration, Fw is the wing aerodynamic force and rCoP is the vector from the fly's center of mass (CoM) to the wing center of pressure (CoP).

rCoP incorporates both the motion of the wing relative to the hinge and a fixed distance between the body CoM and the wing hinge, and can thus be written in the fly body frame as (apostrophe denotes body frame). The vector from the body CoM to the wing hinge, rhinge, is a fixed length estimated from our 3D reconstructions, and is given in the lab frame in Table 2. R is the wing span, and the unit vector corresponding to the direction of the wing span relative to the hinge can be written as:
(A 9)
Note that these are expressions for rCoP, the moment arm vector in a frame with axes fixed to the body; for our calculations, we used a fixed-orientation frame, and thus calculated rCoP=QrCoP, where Q is a rotation matrix. We determine Fw using the quasi-steady model in Eqns 1 and 2 of the main text, with the wingtip velocity written in the lab frame as:
(A 10)

We calculated the aerodynamic force Fw, including the counter-flapping torque (CFT) terms described in Hedrick et al. (2009). Note, however, that our equations of motion do not include drag on the body, which we assumed to be negligible. We chose to solve the non-linear equations of motion, rather than the linearized dynamics in Ristroph et al. (2013).

Control is implemented by adjusting the front stroke angle of the prescribed wing kinematics according to Eqn 3. The Ki and Kp parameters were determined by the fit to the experimental data. The inputs for the controller in each wingbeat – the body pitch angle and pitch velocity – were taken as their mean values during the previous wingbeat. This scheme represents a time delay of 1 wingbeat while avoiding the effects of the inherent small-scale pitch oscillations. Before we apply the perturbation, we let the simulated fly stabilize to a steady-state body pitch angle of 45 deg. The perturbation is then applied, with magnitude roughly corresponding to the accelerations observed in the experiments. In simulation runs that tested controller models other than PI, we let the system stabilize at θb=45 deg using a PI controller and only then applied the perturbation and simultaneously changed the controller type.

Analysis of data from previous studies

Following Beatus et al. (2015), we used the PI control model to predict the pitch response of tethered fruit flies previously published by Dickinson (1999). In Dickinson (1999), flies were tethered to a gimbal apparatus that oscillated about different rotational axes. The left wingstroke amplitude of the flies, ɸleft, was measured using photodetectors that recorded the wings' shadows. As noted by Dickinson (1999), flies do not adjust their back stroke angle during pitch correction, so stroke amplitude is a good proxy for , the quantity we measured in free-flight experiments.

In one of the measurements reported in Dickinson (1999), pitch perturbations were imposed so that: θb(t)=Asin(ωt) and , with A=25 deg, period T=0.63 s and maximum pitch velocity 250 deg s−1. The left wingstroke amplitude was plotted against both pitch angle and pitch velocity (fig. 3A,C in Dickinson, 1999). Using standard image-processing techniques, we extract the data from these plots.

The pitch oscillations in the tethered experiments have period of 630 ms, which is much longer than the observed pitch correction latency times from our experiments (≈10 ms), so we consider the controller time delay negligible. We then write the form for our controller, now in terms of left wingstroke amplitude, as:
(A 11)

where the left wingstroke amplitude, ɸleft(t), and mean stroke amplitude, ɸmean, are related to by the linear relationship: . Considering only the left wing does not reduce the generality of this analysis, as pitch correction is left/right symmetric. We manually fitted for the control parameters from the data. The fitted parameters obtained are Ki=0.3 and Kp=8 ms, comparable to the parameters from the main text.

The predictions of the PI controller fit are shown in Figs S1–3. The output of the PI controller is plotted as a function of both pitch angle and pitch velocity in Fig. S1, yielding an ellipse in both cases (R2=0.761 and 0.842, respectively). The linear model for ɸleft as a function of gives R2=0.556. We also plotted the PI controller prediction in the 3D space whose axes are (θ­b, , ɸleft) (Fig. S2). The PI controller predicts an inclined ellipse in this space, the projections of which onto the horizontal axes yield the plots in Fig. S1. The inclination of the ellipse shows that the corrective response depends on both pitch angle and pitch velocity, i.e. the measured data in Dickinson (1999) are consistent with a PI controller.

Additionally, we show the predicted output of the PI controller plotted as a function of pitch acceleration in Fig. S3. Consistent with Dickinson (1999), Fig. S3 shows that the fly's corrective response can be quantitatively captured without including information about the pitch acceleration. Fig. S3 suggests that pitch acceleration is unimportant in determining the flies' corrective wing kinematics, and thus excludes a PID controller model.

We thank Grace (Li) Chi and Andy Clark for providing flies; Ty Hedrick for advice on camera calibration; and Andy Ruina and the Cohen group members for useful conversations. We especially thank Professor Noah Cowan for his insights and the important role he has played in this and previous manuscripts helping us develop and apply control theory models to our systems.

Author contributions

S.C.W. and L.C. conducted experiments and collected data. T.B. and S.C.W. developed code for data analysis and simulation. S.C.W. analyzed the data. S.C.W., T.B. and I.C. wrote the paper.

Funding

This work was supported in part by a National Science Foundation (NSF) DMR award (no. 1056662) and in part by an Army Research Office (ARO) award (no. 61651-EG). S.C.W. was supported by a National Defense Science and Engineering Graduate (NDSEG) Fellowship. T.B. was supported by the Cross Disciplinary Post-Doctoral Fellowship of the Human Frontier Science Program.

Balint
,
C. N.
and
Dickinson
,
M. H.
(
2004
).
Neuromuscular control of aerodynamic forces and moments in the blowfly, Calliphora vicina
.
J. Exp. Biol.
207
,
3813
-
3838
.
Beatus
,
T. B.
,
Guckenheimer
,
J. M.
and
Cohen
,
I.
(
2015
).
Controlling roll perturbations in fruit flies
.
J. R. Soc. Interface
12
,
20150075
.
Bender
,
J. A.
and
Dickinson
,
M. H.
(
2006
).
Visual stimulation of saccades in magnetically tethered Drosophila
.
J. Exp. Biol.
209
,
3170
-
3182
.
Bergou
,
A. J.
,
Ristroph
,
L.
,
Guckenheimer
,
J.
,
Cohen
,
I.
and
Wang
,
Z. J.
(
2010
).
Fruit flies modulate passive wing pitching to generate in-flight turns
.
Phys. Rev. Lett.
104
,
148101
.
Card
,
G.
and
Dickinson
,
M. H.
(
2008
).
Performance trade-offs in the flight initiation of Drosophila
.
J. Exp. Biol.
211
,
341
-
353
.
Chang
,
S.
and
Wang
,
Z. J.
(
2014
).
Predicting fruit flys sensing rate with insect flight simulations
.
Proc. Natl. Acad. Sci. USA
111
,
11246
-
11251
.
Cheng
,
B.
,
Fry
,
S. N.
,
Huang
,
Q.
,
Dickson
,
W. B.
,
Dickinson
,
M. H.
, and
Deng
,
X.
(
2009
).
Turning dynamics and passive damping in flapping flight
. in
Robotics and Automation, 2009. ICRA ‘09. IEEE International Conference on (IEEE)
pp.
1889
-
1896
.
Cheng
,
B.
,
Deng
,
X.
and
Hedrick
,
T. L.
(
2011
).
The mechanics and control of pitching manoeuvres in a freely flying hawkmoth (Manduca sexta)
.
J. Exp. Biol.
214
,
4092
-
4106
.
Combes
,
S. A.
and
Dudley
,
R.
(
2009
).
Turbulence-driven instabilities limit insect flight performance
.
Proc. Natl. Acad. Sci. USA
106
,
9105
-
9108
.
Cowan
,
N. J.
and
Fortune
,
E. S.
(
2007
).
The critical role of locomotion mechanics in decoding sensory systems
.
J. Neurosci.
27
,
1123
-
1128
.
Cowan
,
N. J.
,
Lee
,
J.
and
Full
,
R. J.
(
2006
).
Task-level control of rapid wall following in the American cockroach
.
J. Exp. Biol.
209
,
1617
-
1629
.
Cowan
,
N. J.
,
Ankarali
,
M. M.
,
Dyhr
,
J. P.
,
Madhav
,
M. S.
,
Roth
,
E.
,
Sefati
,
S.
,
Sponberg
,
S.
,
Stamper
,
S. A.
,
Fortune
,
E. S.
and
Daniel
,
T. L.
(
2014
).
Feedback control as a framework for understanding tradeoffs in biology
.
Integr. Comp. Biol.
54
,
223
-
237
.
Dickerson
,
A. K.
,
Shankles
,
P. G.
,
Madhavan
,
N. M.
and
Hu
,
D. L.
(
2012
).
Mosquitoes survive raindrop collisions by virtue of their low mass
.
Proc. Natl. Acad. Sci. USA
109
,
9822
-
9827
.
Dickinson
,
M. H.
(
1999
).
Haltere-mediated equilibrium reflexes of the fruit fly, Drosophila melanogaster
.
Philos. Trans. R. Soc. B Biol. Sci.
354
,
903
-
916
.
Dickinson
,
M. H.
and
Tu
,
M. S.
(
1997
).
The function of dipteran flight muscle
.
Comp. Biochem. Physiol. A Physiol.
116
,
223
-
238
.
Dickinson
,
M. H.
,
Lehmann
,
F.-O.
and
Sane
,
S. P.
(
1999
).
Wing rotation and the aerodynamic basis of insect flight
.
Science
284
,
1954
-
1960
.
Dyhr
,
J. P.
,
Morgansen
,
K. A.
,
Daniel
,
T. L.
and
Cowan
,
N. J.
(
2013
).
Flexible strategies for flight control: an active role for the abdomen
.
J. Exp. Biol.
216
,
1523
-
1536
.
Elzinga
,
M. J.
,
Dickson
,
W. B.
and
Dickinson
,
M. H.
(
2012
).
The influence of sensory delay on the yaw dynamics of a flapping insect
.
J. R. Soc. Interface
9
,
1685
-
1696
.
Ennos
,
A. R.
(
1989
).
The kinematics and aerodynamics of the free flight of some diptera
.
J. Exp. Biol.
142
,
49
-
85
.
Faruque
,
I.
and
Sean Humbert
,
J.
(
2010
).
Dipteran insect flight dynamics. part 1 longitudinal motion about hover
.
J. Theor. Biol.
264
,
538
-
552
.
Fayyazuddin
,
A.
and
Dickinson
,
M. H.
(
1999
).
Convergent mechanosensory input structures the firing phase of a steering motor neuron in the blowfly, calliphora
.
J. Neurophysiol.
82
,
1916
-
1926
.
Fox
,
J. L.
and
Daniel
,
T. L.
(
2008
).
A neural basis for gyroscopic force measurement in the halteres of holorusia
.
J. Comp. Physiol. A
194
,
887
-
897
.
Fox
,
J. L.
,
Fairhall
,
A. L.
and
Daniel
,
T. L.
(
2010
).
Encoding properties of haltere neurons enable motion feature detection in a biological gyroscope
.
Proc. Natl. Acad. Sci. USA
107
,
3840
-
3845
.
Fry
,
S. N.
,
Sayaman
,
R.
and
Dickinson
,
M. H.
(
2003
).
The aerodynamics of free-flight maneuvers in Drosophila
.
Science
300
,
495
-
498
.
Fry
,
S. N.
,
Sayaman
,
R.
and
Dickinson
,
M. H.
(
2005
).
The aerodynamics of hovering flight in Drosophila
.
J. Exp. Biol.
208
,
2303
-
2318
.
Gao
,
N.
,
Aono
,
H.
and
Liu
,
H.
(
2011
).
Perturbation analysis of 6dof flight dynamics and passive dynamic stability of hovering fruit fly Drosophila melanogaster
.
J. Theor. Biol.
270
,
98
-
111
.
Hafting
,
T.
,
Fyhn
,
M.
,
Molden
,
S.
,
Moser
,
M.-B.
and
Moser
,
E. I.
(
2005
).
Microstructure of a spatial map in the entorhinal cortex
.
Nature
436
,
801
-
806
.
Hedrick
,
T. L.
,
Cheng
,
B.
and
Deng
,
X.
(
2009
).
Wingbeat time and the scaling of passive rotational damping in flapping flight
.
Science
324
,
252
-
255
.
Heide
,
G.
(
1983
).
Neural mechanisms of flight control in Diptera
.
BIONA Rep.
2
,
35
-
52
.
Heide
,
G.
and
Götz
,
K. G.
(
1996
).
Optomotor control of course and altitude in Drosophila melanogaster is correlated with distinct activities of at least three pairs of flight steering muscles
.
J. Exp. Biol.
199
,
1711
-
1726
.
Hesselberg
,
T.
and
Lehmann
,
F.-O.
(
2007
).
Turning behaviour depends on frictional damping in the fruit fly Drosophila
.
J. Exp. Biol.
210
,
4319
-
4334
.
Huston
,
S. J.
and
Krapp
,
H. G.
(
2009
).
Nonlinear integration of visual and haltere inputs in fly neck motor neurons
.
J. Neurosci.
29
,
13097
-
13105
.
Lee
,
J.
,
Sponberg
,
S. N.
,
Loh
,
O. Y.
,
Lamperski
,
A. G.
,
Full
,
R. J.
and
Cowan
,
N. J.
(
2008
).
Templates and anchors for antenna-based wall following in cockroaches and robots
.
IEEE Trans. Robotics
24
,
130
-
143
.
Lin
,
H.-T.
,
Ros
,
I. G.
and
Biewener
,
A. A.
(
2014
).
Through the eyes of a bird: modelling visually guided obstacle flight
.
J. R. Soc. Interface
11
,
20140239
.
Liu
,
H.
,
Nakata
,
T.
,
Gao
,
N.
,
Maeda
,
M.
,
Aono
,
H.
and
Shyy
,
W.
(
2010
).
Micro air vehicle-motivated computational biomechanics in bio-flights: aerodynamics, flight dynamics and maneuvering stability
.
Acta Mech. Sin.
26
,
863
-
879
.
Lourakis
,
M. I. A.
and
Argyros
,
A. A.
(
2009
).
SBA: a software package for generic sparse bundle adjustment
.
ACM Trans. Math. Software
36
,
1
-
30
.
Madhav
,
M. S.
,
Stamper
,
S. A.
,
Fortune
,
E. S.
and
Cowan
,
N. J.
(
2013
).
Closed-loop stabilization of the jamming avoidance response reveals its locally unstable and globally nonlinear dynamics
.
J. Exp. Biol.
216
,
4272
-
4284
.
Miyan
,
J. A.
and
Ewing
,
A. W.
(
1985
).
How diptera move their wings: a re-examination of the wing base articulation and muscle systems concerned with flight
.
Philos. Trans. R. Soc. B Biol. Sci.
311
,
271
-
302
.
Muijres
,
F. T.
,
Elzinga
,
M. J.
,
Melis
,
J. M.
and
Dickinson
,
M. H.
(
2014
).
Flies evade looming targets by executing rapid visually directed banked turns
.
Science
344
,
172
-
177
.
Nalbach
,
G.
(
1994
).
Extremely non-orthogonal axes in a sense organ for rotation: behavioural analysis of the dipteran haltere system
.
Neuroscience
61
,
149
-
163
.
Ofstad
,
T. A.
,
Zuker
,
C. S.
and
Reiser
,
M. B.
(
2011
).
Visual place learning in Drosophila melanogaster
.
Nature
474
,
204
-
207
.
Ortega-Jimenez
,
V. M.
,
Greeter
,
J. S. M.
,
Mittal
,
R.
and
Hedrick
,
T. L.
(
2013
).
Hawkmoth flight stability in turbulent vortex streets
.
J. Exp. Biol.
216
,
4567
-
4579
.
Pérez-Arancibia
,
N. O.
,
Ma
,
K. Y.
,
Galloway
,
K. C.
,
Greenberg
,
J. D.
and
Wood
,
R. J.
(
2011
).
First controlled vertical flight of a biologically inspired microrobot
.
Bioinspir. Biomim.
6
,
036009
.
Pringle
,
J. W. S.
(
1948
).
The gyroscopic mechanism of the halteres of Diptera
.
Philos. Trans. R. Soc. B Biol. Sci.
233
,
347
-
384
.
Pringle
,
J. W. S.
(
1949
).
The excitation and contraction of the flight muscles of insects
.
J. Physiol.
108
,
226
-
232
.
Ravi
,
S.
,
Crall
,
J. D.
,
Fisher
,
A.
and
Combes
,
S. A.
(
2013
).
Rolling with the flow: bumblebees flying in unsteady wakes
.
J. Exp. Biol.
216
,
4299
-
4309
.
Ristroph
,
L.
,
Berman
,
G. J.
,
Bergou
,
A. J.
,
Wang
,
Z. J.
and
Cohen
,
I.
(
2009
).
Automated hull reconstruction motion tracking (hrmt) applied to sideways maneuvers of free-flying insects
.
J. Exp. Biol.
212
,
1324
-
1335
.
Ristroph
,
L.
,
Bergou
,
A. J.
,
Ristroph
,
G.
,
Coumes
,
K.
,
Berman
,
G. J.
,
Guckenheimer
,
J.
,
Wang
,
Z. J.
and
Cohen
,
I.
(
2010
).
Discovering the flight autostabilizer of fruit flies by inducing aerial stumbles
.
Proc. Natl. Acad. Sci. USA
107
,
4820
-
4824
.
Ristroph
,
L.
,
Bergou
,
A. J.
,
Guckenheimer
,
J.
,
Wang
,
Z. J.
and
Cohen
,
I.
(
2011
).
Paddling mode of forward flight in insects
.
Phys. Rev. Lett.
106
,
178103
.
Ristroph
,
L.
,
Ristroph
,
G.
,
Morozova
,
S.
,
Bergou
,
A. J.
,
Chang
,
S.
,
Guckenheimer
,
J.
,
Wang
,
Z. J.
and
Cohen
,
I.
(
2013
).
Active and passive stabilization of body pitch in insect flight
.
J. R. Soc. Interface
10
,
20130237
.
Roth
,
E.
,
Zhuang
,
K.
,
Stamper
,
S. A.
,
Fortune
,
E. S.
and
Cowan
,
N. J.
(
2011
).
Stimulus predictability mediates a switch in locomotor smooth pursuit performance for Eigenmannia Virescens
.
J. Exp. Biol.
214
,
1170
-
1180
.
Roth
,
E.
,
Sponberg
,
S.
and
Cowan
,
N. J.
(
2014
).
A comparative approach to closed-loop computation
.
Curr. Opin. Neurobiol.
25
,
54
-
62
.
Sane
,
S. P.
and
Dickinson
,
M. H.
(
2001
).
The control of flight force by a flapping wing: lift and drag production
.
J. Exp. Biol.
204
,
2607
-
2626
.
Seelig
,
J. D.
and
Jayaraman
,
V.
(
2013
).
Feature detection and orientation tuning in the Drosophila central complex
.
Nature
503
,
262
-
266
.
Sherman
,
A.
and
Dickinson
,
M. H.
(
2003
).
A comparison of visual and haltere-mediated equilibrium reflexes in the fruit fly Drosophila melanogaster
.
J. Exp. Biol.
206
,
295
-
302
.
Sherman
,
A.
and
Dickinson
,
M. H.
(
2004
).
Summation of visual and mechanosensory feedback in Drosophila flight control
.
J. Exp. Biol.
207
,
133
-
142
.
Sun
,
M.
(
2014
).
Insect flight dynamics: stability and control
.
Rev. Mod. Phys.
86
,
615
-
646
.
Sun
,
M.
and
Xiong
,
Y.
(
2005
).
Dynamic flight stability of a hovering bumblebee
.
J. Exp. Biol.
208
,
447
-
459
.
Sun
,
M.
,
Wang
,
J.
and
Xiong
,
Y.
(
2007
).
Dynamic flight stability of hovering insects
.
Acta Mech. Sin.
23
,
231
-
246
.
Tammero
,
L. F.
and
Dickinson
,
M. H.
(
2002
).
The influence of visual landscape on the free flight behavior of the fruit fly Drosophila melanogaster
.
J. Exp. Biol.
205
,
327
-
343
.
Taylor
,
G. K.
(
2001
).
Mechanics and aerodynamics of insect flight control
.
Biol. Rev.
76
,
449
-
471
.
Taylor
,
G. K.
and
Thomas
,
A. L. R.
(
2003
).
Dynamic flight stability in the desert locust Schistocerca gregaria
.
J. Exp. Biol.
206
,
2803
-
2829
.
Taylor
,
G. K.
and
Żbikowski
,
R. I.
(
2005
).
Nonlinear time-periodic models of the longitudinal flight dynamics of desert locusts Schistocerca gregaria
.
J. R. Soc. Interface
2
,
197
-
221
.
Theriault
,
D. H.
,
Fuller
,
N. W.
,
Jackson
,
B. E.
,
Bluhm
,
E.
,
Evangelista
,
D.
,
Wu
,
Z.
,
Betke
,
M.
and
Hedrick
,
T. L.
(
2014
).
A protocol and calibration method for accurate multi-camera field videography
.
J. Exp. Biol.
217
,
1843
-
1848
.
Tytell
,
E. D.
,
Holmes
,
P.
and
Cohen
,
A. H.
(
2011
).
Spikes alone do not behavior make: why neuroscience needs biomechanics
.
Curr. Opin. Neurobiol.
21
,
816
-
822
.
Vance
,
J. T.
,
Faruque
,
I.
and
Humbert
,
J. S.
(
2013
).
Kinematic strategies for mitigating gust perturbations in insects
.
Bioinspir. Biomim.
8
,
016004
.
Walker
,
S. M.
,
Schwyn
,
D. A.
,
Mokso
,
R.
,
Wicklein
,
M.
,
Müller
,
T.
,
Doube
,
M.
,
Stampanoni
,
M.
,
Krapp
,
H. G.
and
Taylor
,
G. K.
(
2014
).
In vivo time-resolved microtomography reveals the mechanics of the blowfly flight motor
.
PLoS Biol.
12
,
e1001823
.
Warrick
,
D.
,
Hedrick
,
T.
,
Fernández
,
M. J.
,
Tobalske
,
B.
and
Biewener
,
A.
(
2012
).
Hummingbird flight
.
Curr. Biol.
22
,
R472
-
R477
.
Windsor
,
S. P.
,
Bomphrey
,
R. J.
and
Taylor
,
G. K.
(
2014
).
Vision-based flight control in the hawkmoth
.
Hyles lineata. J. R. Soc. Interface
11
,
20130921
.
Xu
,
N.
and
Sun
,
M.
(
2013
).
Lateral dynamic flight stability of a model bumblebee in hovering and forward flight
.
J. Theor. Biol.
319
,
102
-
115
.
Yartsev
,
M. M.
,
Witter
,
M. P.
and
Ulanovsky
,
N.
(
2011
).
Grid cells without theta oscillations in the entorhinal cortex of bats
.
Nature
479
,
103
-
107
.
Zanker
,
J. M.
(
1988
).
On the mechanism of speed and altitude control in drosophila melanogaster
.
Physiol. Entomol.
13
,
351
-
361
.
Zanker
,
J. M.
(
1990
).
The wing beat of Drosophila melanogaster. III. control
.
Philos. Trans. R. Soc. B Biol. Sci.
327
,
43
-
64
.
Zhang
,
Y.
and
Sun
,
M.
(
2010
).
Dynamic flight stability of a hovering model insect: lateral motion
.
Acta Mech. Sin.
26
,
175
-
190
.
Zhang
,
Y.-L.
and
Sun
,
M.
(
2011
).
Stabilization control of a hovering model insect: lateral motion
.
Acta Mech. Sin.
27
,
823
-
832
.

Competing interests

The authors declare no competing or financial interests.

Supplementary information