Analysis of June 2, 2016 bolide event
Csaba Palotai, Ramanakumar Sankar, Dwayne L. Free, J. Andreas Howell, Elena Botella, Daniel Batcheldor
MMNRAS , 1–13 (2017) Preprint 22 May 2019 Compiled using MNRAS L A TEX style file v3.0
Analysis of the June 2, 2016 bolide event over Arizona
Csaba Palotai, (cid:63) Ramanakumar Sankar, Dwayne L. Free, J. Andreas Howell, Elena Botella and Daniel Batcheldor Department of Physics & Space Sciences, Florida Institute of Technology, 150 W. University Blvd., Melbourne, FL 32901, USA Spalding Allsky Camera Network, SkySentinel, LLC, 958 Shaw Circle, Melbourne, FL 32940, USA Department of Online Science, Florida Institute of Technology, 150 W. University Blvd., Melbourne, FL 32901, USA
Accepted XXX. Received YYY; in original form ZZZ
ABSTRACT
On June 2, 2016 at 10h56m UTC, a − . ± . magnitude superbolide was observedover Arizona. Fragments were located a few days later and the meteorites were giventhe name Dishchii’bikoh. We present analysis of this event based on 3 cameras and amulti-spectral sensor observations by the SkySentinel continuous fireball-monitoringcamera network, supplemented by a dash cam footage and a fragmentation model.The bolide began its luminous flight at an altitude of . ± . km at coordinates φ = . ± . ° N planetographic latitude and λ = . ± . ° W longitude, andit had a pre-atmospheric velocity of . ± . km/s. The calculated orbital parametersindicate that the meteoroid did not belong to any presently known asteroid family.From our calculations, the impacting object had an initial mass of . ± . metrictonnes with an estimated initial diameter of . ± . m. Key words: meteors – meteor light curve – asteroids
At 10h 56m 27s UTC on June 2, 2016 a bright fireball wasobserved over Arizona. The American Meteor Society re-ceived 421 reports about this fireball, most of them from Ari-zona but people also witnessed this event in Utah, New Mex-ico, California, Texas, Colorado and Nevada. Sonic boomsfrom the bolide were heard across the greater Phoenix area.Videos of the bolide from dash cams and security camerasappeared on YouTube and various media outlets. After sun-rise, videos captured the dust trail that the impactor leftbehind.Camera and satellite observations of meteors have beenused for decades in order to help to determine the mass,trajectory and orbital parameters of the impacting body.Various satellites, NASA’s All Sky Fireball Network andthe Lowell Observatory Cameras for All-Sky Meteor Surveil-lance (LO-CAMS) also recorded the Arizona bolide event.The Arizona Geological Survey’s seismic network picked upa signal near Payson that was consistent with an airburstevent and according to the agency it marked the explosionof the asteroid. This was the largest observed bolide eventover the continental United States since March 2010 thatwas reported by the Center for Near Earth Object Studies (cid:63)
E-mail: cpalotai@fit.edu (CsP) (CNEOS) , and the fifth largest since 1988 when the agencybegan recording fireball data from US Government sensors.The bolide was also observed by multiple nodes of theSkySentinel’s Spalding Allsky Camera Network (Figure 1)and the recorded data allow us to characterize some of thephysical properties and the likely origin of the impactingbody. Here, we present our analysis of the available SkySen-tinel observations, coupled with complementary data andmodeling tools. We also compare our inferred results withdata reported by other sources. Mr. R.E. Spalding (1936-2017, Sandia National Labs) de-veloped the Allsky Camera System to monitor, track, andanalyze large meteor events in order to provide “ground-truth” to assist both science (NASA) and treaty monitor-ing (Nuclear-Test-Ban Treaty Organization - CTBTO) op-erations in confirming the impact of large meteor fireballs(bolides) in Earth’s atmosphere. The collected data is alsoused to support the refinement of the energy calculations of h ttps://cneos.jpl.nasa.gov/fireballs/, accessed December 2,2017 © a r X i v : . [ a s t r o - ph . E P ] M a y Cs. Palotai et al. those events, as well as, to improve trajectory calculationsand orbit determination for the impacting bodies. Anotherimportant goal of the project was to develop a companioninstrument, MultiSpectral Radiometer (MSR) for compari-son to the Allsky cameras, and government-collected data,in order to improve the diagnostic capability of SkySentinel.The network consists of a large number of wide-angleview cameras at various sites throughout the continentalUnited States and other countries and it also includes theinfrastructure that permits the archival of the observationaldata that is available for processing and analysis for the sci-entific community. Software tools were developed and addedfor calibration, removal of detector effects and anomalies, au-tomatic event detection and correlation among stations, andautomatic trajectory computation.Mr. R. E. Spalding formally established the current Sky-Sentinel Allsky Camera Network that was transitioned toSkySentinel, LLC and is now operated as a Joint FloridaInstitute of Technology and SkySentinel, LLC Science Edu-cation and Research Program. This paper is a demonstrationof improvements made to SkySentinel under the Joint Pro-gram. The current system was renamed the Spalding AllskyCamera Network in honor of its founder. Data for recordedevents and details on the camera system can be found at . The 2 June 2016 Arizona bolide event was detected byseven SkySentinel nodes. Table 1 lists the SkySentinel nodesthat observed the event. These nodes include a Sony HB-710E Starlight B/W CCD camera coupled with a FujinonYV22X14ASA 3.4mm lens. The cameras are sensitive downto about + visual magnitude using 8.5 second image com-posites. For individual frames recorded at 30 frames per sec-ond, the cameras are sensitive to just above magnitude .The node cameras are enabled throughout the nightand the video feed goes directly into the WSentinel trackingsoftware which automatically triggers on motion/thresholdcrossing and saves 1s of presamples and the entirety ofthe event in four files: composite .jpg image, .mp4 video(640x480 @ 30fps; 8 bits per pixel), .csv spreadsheet of timestamped azimuth/elevation, sum of pixels triggering, sum ofenergy above threshold per pixel, X/Y pixel map of the eventand a .txt file of the start/stop time of the event recording.WSentinel also takes an hourly composite .png image (ad-justable, 8.5s minimum) that is used to calibrate the cameraby star mapping tie points. SkySentinel nodes report positional data that are expressedin Cartesian coordinates ( x , y ) . During the calibration of thedata we convert each SkySentinel camera’s reported ( x , y ) measurements into horizon coordinates ( z , a ) , where z is thezenith distance and a is the azimuth. The image calibra-tion is accomplished by fitting a 10-parameter mathemati-cal model to observations of calibration stars in each node’sCCD imagery. Table 1 lists the SkySentinel nodes that werecalibrated, the dates of the SkySentinel images used for cal- Figure 1.
Composite images of the bolide from Node 6 (Flagstaff,AZ) and Node 79 (Turkey Springs, AZ). Only the first 3s wereused from each video to create the composites due to saturationof subsequent frames. The event itself continues for about an ad-ditional s and brightens by several orders of magnitude. ibration, and the number of calibration points (stars) em-ployed.The method used to perform the calibration was basedon a geometric solution originally developed by Borovickaet al. (1995) to calibrate emulsion-based plates from all-skycameras. Boroviˇcka’s complete solution, which employed 13calibration parameters, achieves a residual error of just 0.015degree. The solution’s main drawback is that it requires a“large number” of calibration stars, including stars at largezenith angles, which are “seldom available”.The Boroviˇcka solution uses a coordinate system thatdiffers in an important way from SkySentinel’s coordinatesystem. To distinguish between them, SkySentinel coordi-nates will henceforth be denoted using lower-case letters ( x , y ) , and upper-case letters will denote Boroviˇcka coordi-nates ( X , Y ) . The difference is that coordinates are inter-changed, i.e., ( X , Y ) correlates to ( y , x ) . For example, ( x = , y = ) in SkySentinel coordinates is the same point as ( X = , Y = ) in Boroviˇcka coordinates.On the other hand, the two coordinate systems are alikein that they both define the center of the upper left pixel tobe the origin of the coordinate system ( , ) . For SkySentinelcameras aligned approximately with North up, the origin ( , ) is the center of the pixel in the Northeast corner of theimage.New Mexico State University, through grant funding,studied how to accomplish calibration of CCD-based videometeor cameras used in SkySentinel (Bannister et al. 2013).Bannister made several simplifications to the Boroviˇckamodel, one of which was to replace the exponential modelof zenith distance by a quadratic approximation. Appliedto actual SkySentinel sensors, Bannister showed that his 8-parameter all-sky calibration model could achieve 0.21 de-gree accuracy in azimuth and 0.09 degree in zenith angle.For a SkySentinel sensor, this is sub-pixel precision because1 pixel subtends about 0.3 degree on the sky.The calibration of SkySentinel sensors is based uponBannister’s paper, but with two significant changes. First,an exponential model replaced the quadratic model of zenithdistance. This is because the quadratic model of zenith dis-tance was leading to systematic error at zenith distancesgreater than 70 ° . Most bolides are seen at large zenith dis-tances (i.e., low elevations), and utmost accuracy is neededin this regime. Indeed, replacement of the quadratic model MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Table 1.
Planetographic coordinates for the SkySentinel nodes that observed the event, and the dates and number of stars used for theastrometric calibration of the pixel positions from the cameras. Nodes 5, 6, 8, 37 and 79 are used in this study.Node Location Imagery Dates No. Cal Stars Latitude ( ° ) Longitude ( ° ) Altitude (m)1 Las Cruces, NM Jun 1-4 130 32.281 -106.754 11912 Las Cruces, NM May 31-Jun 2 54 32.281 -106.754 11915 Alberquerque, NM May 29 - Jun 18 102 35.152 -106.556 16006 Flagstaff, AZ May 31-Jun 2 124 35.200 -111.655 210614 Los Alamos, NM May 31-Jun 4 56 35.887 -106.277 217737 Parker, AZ May 31-Jun 2 109 34.144 -114.291 12879 Payson, AZ May 31-Jun 3 120 34.233 -111.301 15148 ∗ Alberquerque, NM - - 34.951 -106.460 1958 ∗ MultiSpectral Radiometer (MSR) by an exponential model of zenith distance was a stepback to the original Boroviˇcka model (Boroviˇcka (1992);Borovicka et al. (1995)).Second, we discovered that the USB video capture de-vices used to digitize SkySentinel’s analog signal often intro-duced elliptical distortion. Figure 2 illustrates the differencebetween a normal image and one that has elliptical distor-tion. The axes of the ellipse always appeared to align withthe x and y axes of the sensor. The distortion was leading,in some cases, to several degrees periodic error in the resid-uals. To restore image circularity, x and y coordinates wereadjusted using appropriate scale factors. Nodes 5 and 6 hadmeasurable amounts of image ellipticity while Nodes 37 and79 had virtually none. Table 2 tabulates the parameters usedto calibrate each node. Calibration points (stars) brighter than visual magnitude+2 were selected from a list downloaded from the U.S. NavalObservatory’s Online Astronomical Almanac . Hourly Sky-Sentinel composite images were selected from a period of4-5 days centered on the 2 June 2016 bolide event. Imageswere converted to TIF format, 8-bit depth, and enhancedby subtraction of a median combined image to enhance starvisibility. Calibration stars were selected over a wide rangeof azimuths and zenith distances, limited only by availabil-ity of stars bright enough to be seen by SkySentinel. The15-degree/hour change of hour angle made it possible tomeasure the same star multiple times during the course of anight.Mira Pro x64 Windows astronomical software (Mira-metrics Inc. 2017) was used to measure calibration stars un-til there were 50-150 data points. Mira Pro defines the upperleft (Northeast) pixel to have coordinates (1,1), which dif-fers from SkySentinel’s convention that this same pixel hascoordinates (0,0). Accordingly, Mira Pro coordinates weredecremented by 1, along each axis, before data entry intothe SkySentinel calibration model.A Microsoft Excel 2011 workbook was created to solvefor the 10 calibration parameters using a least-squares pro-cedure. For each calibration point, the great circle error be- http://asa.usno.navy.mil/SecH/BrightStarsSearch.html ,accessed July 01, 2015 Figure 2.
Image flattening. The left frame shows a flattening of Φ = . , while the right frame shows Φ = . tween cataloged (USNO Bright Stars) and calculated (Cal-ibration Model) star coordinates was calculated. The errorswere squared and summed to yield a sum of squared er-rors (SSE). Excel’s non-linear solver was then used to min-imize SSE by varying the 10 calibration parameters, untilthe solver converged on a solution. Of these 10 parameters, 9were independent, because the y scale factor was constrainedto equal the reciprocal of the x scale factor. This constraintensured that the area of the image projection remained con-stant. Table 2 lists the values of the calibration parametersfound for each of the three SkySentinel nodes that were usedfor the analysis of the event.Standard deviation of the great circle error, σ = (cid:112) SSE /( n − ) , where n is the number of calibration stars,ranged from 0.07-0.16 degree. Goodness-of-fit was verified byexamining residual plots to ensure they were randomly dis-tributed with no patterns or trends. The mean absolute de-viation of zenith distance | δ z | ranged from 0.04-0.09 degree,and the mean absolute deviation of azimuth | δ a sin z | rangedfrom 0.03-0.07 degree. These uncertainties, no greater thanone-third of a pixel, show that SkySentinel cameras are ca-pable of accurate calibration.In summary, the 10-parameter calibration model pro-duces accurate horizon coordinates ( z , a ) where z is the zenithdistance and a is the azimuth measured from cardinal South.To get azimuth measured from cardinal North, compute a + π , and subtract multiples of 2 π , as needed, to get a resulton the interval [ , π ) . MNRAS000
Image flattening. The left frame shows a flattening of Φ = . , while the right frame shows Φ = . tween cataloged (USNO Bright Stars) and calculated (Cal-ibration Model) star coordinates was calculated. The errorswere squared and summed to yield a sum of squared er-rors (SSE). Excel’s non-linear solver was then used to min-imize SSE by varying the 10 calibration parameters, untilthe solver converged on a solution. Of these 10 parameters, 9were independent, because the y scale factor was constrainedto equal the reciprocal of the x scale factor. This constraintensured that the area of the image projection remained con-stant. Table 2 lists the values of the calibration parametersfound for each of the three SkySentinel nodes that were usedfor the analysis of the event.Standard deviation of the great circle error, σ = (cid:112) SSE /( n − ) , where n is the number of calibration stars,ranged from 0.07-0.16 degree. Goodness-of-fit was verified byexamining residual plots to ensure they were randomly dis-tributed with no patterns or trends. The mean absolute de-viation of zenith distance | δ z | ranged from 0.04-0.09 degree,and the mean absolute deviation of azimuth | δ a sin z | rangedfrom 0.03-0.07 degree. These uncertainties, no greater thanone-third of a pixel, show that SkySentinel cameras are ca-pable of accurate calibration.In summary, the 10-parameter calibration model pro-duces accurate horizon coordinates ( z , a ) where z is the zenithdistance and a is the azimuth measured from cardinal South.To get azimuth measured from cardinal North, compute a + π , and subtract multiples of 2 π , as needed, to get a resulton the interval [ , π ) . MNRAS000 , 1–13 (2017)
Cs. Palotai et al.
Table 2.
Calibration Parameters of three SkySentinel Nodes thatwere used for analysis SkySentinel Node5 6 37 79V 0.005891 0.00589 0.005792 0.005934S 0.002081 0.001458 0.001791 0.001979D 0.017449 0.0191 0.017863 0.017886a [deg] -22.349 -1.14 -2.116 -17.739E [deg] 240.671 59.909 223.085 289.919 (cid:15) [deg] 2.799 2.588 1.449 0.954X x scale factor 0.988604 0.988506 1.000282 0.999640 y scale factor 1.011527 1.011627 0.999718 1.000360 For the trajectory analysis of the bolide we used data fromthree sites: Flagstaff (Node 6), Payson (Node 79), both lo-cated in Arizona, and Albuquerque (Node 5), in New Mex-ico. We could not use data from the other cameras (Nodes1, 2 and 14), as the bolide disappeared below the horizonor behind obstructing objects causing only the bloom to beseen in the videos. In the case of Node 37, there was anunknown debris covering the camera near the location ofthe bolide. This led to an inaccurate centroid for the bolideon the frame, which skewed the results for the trajectoryfrom this camera. By carrying out sensitivity tests, later wefound that increasing the calculated elevation of the bolideby about . ° fixed the offset, but due to the nature of theissue, it is not possible to precisely determine this value. Assuch, we have also excluded Node 37 data from the analysis.The Flagstaff station was the first to detect light fromthis event and in this paper we use the first detection ofthe bolide from that video as a reference time. The Multi-Spectral Radiometer (MSR) recorded data at 10,000 framesper second and used GPS clocks for timing. Consequently,we used the MSR to calibrate the timing in the other ob-servations by aligning the features from the respective lightcurves. The offsets in the Node 79, Node 6 and Node 37data were t = . s, t = . s and t = . s respectively.Hereafter, all times are quoted with reference to the start ofthe event at 10h 56m 27.489s UTC.We used frames from the first few seconds of the record-ings up until the saturation of the pixels made it impossibleto accurately determine the location of the bolide from thecamera frames. This gave us 74 frames from Nodes 6 and57 frames from Node 79. For Node 5, we were able to useframes after t ∼ . s – i.e. 46 frames. Therefore, the first . sof the trajectory was determined only using Nodes 6 and 79,while the later part was determined using all three cameras,until Node 79 saturates. The last . s of the triangulatedpart of the trajectory are determined using Nodes 5 and 6. In order to reconstruct the trajectory we turned to theplane intersection method detailed in Ceplecha (1987). Dueto the locations of our cameras relative to the direction ofmotion of the object, the angular separation between the in-tersecting planes from two of the cameras (Nodes 6 and 79)was small enough ( ∼ . ° ) to result in loss of their statisticalsignificance. This made the plane intersection method im-practical for the analysis, as the statistical weight of this pairis only about . For all three cameras, however, we wereable to use an alternate approach that employs the least-squares technique, similar to that described by Borovicka(1990).Using this triangulation method, for each time step wedefine a vector that points to the location of the bolide in ageocentric reference frame. Here, R i is the unit vector point-ing to the observation site from the center of the Earth and r i is the unit vector pointing to the bolide from the observa-tion site. The location of the bolide is then the vector sumwith R i being the distance from the observation site to thecenter of the Earth and r i being the free parameter whichis the distance between the observation site and the bolide.The location of the bolide in terms of the site i is then: h i = r i r i + R i R i (1)The solution to the triangulation is obtained by finding r i and r j which minimizes the distance between the positionof the bolide from any pair of observations i , j , using theGauss-Newton method. We define the difference vector d : d ( r i , r j ) = h i − h j = r i r i + R i R i − r j r j − R j R j (2)Taking u , v and w to be the x-, y- and z- coordinates of d , we can create the Jacobian J and b , given by. J = dudr i dudr j dvdr i dvdr j dwdr i dwdr j b = − u ( r i , r j ) v ( r i , r j ) w ( r i , r j ) (3)Choosing an initial guess for p = ( r i , r j ) , the solu-tion is found iteratively by solving for the step given by d p = ( J T J ) − J T b , until the | d | is below a threshold, or | d p | = . The pairs of observation each produce a solutionto the meteor position as a function of time, i.e. in our case,for 3 observations, we have 6 solutions. A standard linearregression is applied along each direction for all solutions,and the corresponding residuals are determined.The details of the trajectory parameters are listed in Ta-ble 3. The luminous flight began at coordinates φ = . ± . ° N planetographic latitude, λ = . ± . ° W lon-gitude, and a height of . ± . km.To obtain the radiant, we corrected the velocity forEarth’s rotation and the curvature of the trajectory due toEarth’s gravity, following the method described in Ceplecha(1987). This was then projected onto the celestial spheredefined by the geocentric inertial frame. The geocentric ve-locity was determined at the entry location using v g = v ∞ − GM ⊕ R ⊕ + h b , (4)where M ⊕ and R ⊕ is the mass of the Earth and the WGS84radius at the current location respectively, and h b is the MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Node 79Node 6 Node 5Node 37 MSRNode 1, 2Node 14TucsonPhoenix Albuquerque
AZ NM
Figure 3.
The trajectory of the bolide with the entry locationmarked with a black circle. Coloured circles indicate the camerasthat we used to determine the trajectory. The stars represent thecameras that were only used for light curve analysis, while theblack crosses denote cameras that observed the bolide but werenot used in the calculations. Nearby cities are marked with thesmall black circles. Locations of the fragments found are markedby the red × ’s. The gray circle is the location of the car whosedash cam video was used in the analysis. The arrow indicates thedirection that the car was travelling. The black diamond showsthe peak brightness location from CNEOS database. height of the bolide at the entry location. v ∞ is the pre-atmospheric entry speed in the geocentric inertial referenceframe.We also obtained the AZ/EL of the bolide radiant in theco-rotating geocentric frame at the entry point of the bolideusing the AstroPy module in Python (Astropy Collaborationet al. 2018). The Arizona fireball moved closely from Northto South, with the azimuth of the radiant being ± ° ,measured from North.From the SkySentinel observations we could not deter-mine the altitude of the maximum brightness or the terminalheight of the luminous flight due to image saturation and theobscured view that we described above.Figure 3 illustrates the trajectory, along with a mapof the nodes which observed the event. Figure 4 shows thedistance travelled by the bolide along the track as a functionof time while Figure 5 and 6 show the scatter in the results.Node 5 has the largest scatter, although this is likely dueto it being much further away compared to the other nodes.The average deviation is about m. Performing a linear fit to positions retrieved by triangulat-ing SkySentinel camera data between 0.4s and 2.5s gives anobserved velocity of . ± . km/s with respect to the ro-tating surface of the Earth. During this time interval, thereis no observed deceleration of the bolide. Thus, we considerthis value as the pre-atmospheric entry velocity that remainsconstant until the initial breakup.We could not include data for the later stages; duringthe main flare event the frames were saturated and after de-saturation, the bolide disappeared behind surface objects,e.g., a tree in the case of Node 79. This made it impossible Time [s] D i s t a n ce t r a v e ll e d a l o n g t r a c k [ k m ] Figure 4.
The distance travelled by the bolide as a function oftime. The red points near the start of the event are from the tri-angulation of data SkySentinel Nodes 6 and 79, the green crosses( + ) are using all three cameras, the dark-red × ’s are using Nodes5 and 6, while the blue points are from the dash cam video. Thesolid line is from the fragmentation model. to determine from SkySentinel data, the change in velocityduring and after the breakup of the object. For the lightcurve analysis, we attempted photometric cal-ibration of the SkySentinel cameras. These datasets werefrom 8-bit grayscale cameras in the visible wavelength, ob-serving at 30 frames per second. The bolide was below thehorizonal throughout the entire video for Nodes 1, 2, and 14,while for Node 5, the bolide was visible only for the first twoseconds. Therefore, only three cameras (Nodes 6, 37 and 79)observed the entire event.The sum of the raw pixel counts were taken to be pro-portional to the intensity of the light observed, and thus tak-ing the base-10 logarithm and multiplying by a factor of -2.5provided an uncalibrated instrumental magnitude. A photo-metric calibration of the SkySentinel was attempted usingthe Moon and Vega. An image of the Moon was used on alater date for all four cameras. However, the Moon’s lightsaturated the camera, making accurate calibration difficult.In the case of nodes 6 and 37, Vega was visible, but since theSNR of Vega was ∼ in both cameras, we were unable toreliably determine the photometric correction. Furthermore,the raw signal from three cameras (Nodes 6, 37 and 79) areclearly saturated (Fig 7), so the SkySentinel cameras werenot used in the photometry analysis.We also used data from the MSR instrument in NewMexico, which observed at 5 bands: UV, Blue, Green, Redand Infrared at 10,000 frames per second. The quantum ef-ficiency and responsivity of the MSR were already known(Figure 8) and thus, it was a straightforward matter of con-verting the raw signal from the MSR into a calibrated mag- MNRAS000
The distance travelled by the bolide as a function oftime. The red points near the start of the event are from the tri-angulation of data SkySentinel Nodes 6 and 79, the green crosses( + ) are using all three cameras, the dark-red × ’s are using Nodes5 and 6, while the blue points are from the dash cam video. Thesolid line is from the fragmentation model. to determine from SkySentinel data, the change in velocityduring and after the breakup of the object. For the lightcurve analysis, we attempted photometric cal-ibration of the SkySentinel cameras. These datasets werefrom 8-bit grayscale cameras in the visible wavelength, ob-serving at 30 frames per second. The bolide was below thehorizonal throughout the entire video for Nodes 1, 2, and 14,while for Node 5, the bolide was visible only for the first twoseconds. Therefore, only three cameras (Nodes 6, 37 and 79)observed the entire event.The sum of the raw pixel counts were taken to be pro-portional to the intensity of the light observed, and thus tak-ing the base-10 logarithm and multiplying by a factor of -2.5provided an uncalibrated instrumental magnitude. A photo-metric calibration of the SkySentinel was attempted usingthe Moon and Vega. An image of the Moon was used on alater date for all four cameras. However, the Moon’s lightsaturated the camera, making accurate calibration difficult.In the case of nodes 6 and 37, Vega was visible, but since theSNR of Vega was ∼ in both cameras, we were unable toreliably determine the photometric correction. Furthermore,the raw signal from three cameras (Nodes 6, 37 and 79) areclearly saturated (Fig 7), so the SkySentinel cameras werenot used in the photometry analysis.We also used data from the MSR instrument in NewMexico, which observed at 5 bands: UV, Blue, Green, Redand Infrared at 10,000 frames per second. The quantum ef-ficiency and responsivity of the MSR were already known(Figure 8) and thus, it was a straightforward matter of con-verting the raw signal from the MSR into a calibrated mag- MNRAS000 , 1–13 (2017)
Cs. Palotai et al.
Table 3.
Trajectory and orbital parameters using triangulation of the bolide’s position from the initial frames of SkySentinel cameras,the analysis of the dash cam video and the fragmentation model. Azimuth and zenith distance are given for the apparent radiant in theco-rotating geocentric frame. The location of the geocentric radiant is given on the celestial sphere with corrections applied for rotationof the Earth and curvature of the trajectory by Earth’s gravity.Trajectory OrbitEntry Latitude ∗ ± ° N Semi-major axis 1.13 ± ∗ ± ° W Eccentricity 0.210 ± ∗ ± ± ∗ ± ± ° Peak brightness magnitude † -20.4 ± ± ° Height of peak brightness † ± ± ° End height † ± † . ° N Apparent radiant at entry locationEnd longitude † . ° W Zenith distance of radiant ∗ ± ° Azimuth of radiant ± ° Geocentric radiant (J2016.5)Right ascension ± ° Declination ± ° Geocentric velocity . ± . km/s ∗ from triangulation using SkySentinel cameras † using data from dash cam video and the fragmentation model. . . . . . − . − . − . . . . L e n g t h r e s i du a l [ k m ] . . . . . Time [s] − . − . − . . . . L a t e r a l d e v i a t i o n [ k m ] Figure 5.
The residual in the expected length along track (top),and the lateral deviation from the average trajectory (bottom)from the triangulation solution using SkySentinel cameras. Node6 are plotted with × ’s, Node 79 with points and Node 5 withsquares. The colors represent different pairs used for the solution:black for Nodes 6 and 79, red for Nodes 6 and 5, and blue forNodes 79 and 5. . . . . . − . − . − . . . . D e v i a t i o n i n x [ k m ] . . . . . Time [s] − . − . − . . . . D e v i a t i o n i n z [ k m ] . . . . . − . − . − . . . . D e v i a t i o n i n y [ k m ] Figure 6.
The residuals along the x , y and z -direction from thelinear fit. The color coding and markers follow the scheme givenin Fig 5. nitude, from which we obtained the luminosity. We take apower of W at km to correspond to an absolutemagnitude of of zero (Ceplecha et al. 1998). The magnitudevalues were corrected for extinction at each time stamp, withthe extinction coefficients determined from calibrating an
MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Time [s] S i g n a l Node 79Node 6Node 37
Figure 7.
The signal from the SkySentinel cameras in arbitraryunits. All three cameras were saturated during the main flareevent, making the data unusable.
300 400 500 600 700 800 900 1000 1100 1200
Wavelength (nm) . . . . . R e s p o s i v i t y QEUVBlueGreenRedIR
Figure 8.
The dashed curve shows the quantum efficiency of thedetector, while the solid lines are the responsivity (A/W) of eachfilter - from left to right: UV, blue, green, red, IR image of the Moon at different airmass. Figure 9 plots theabsolute magnitudes and luminosity from the MSR data.After the initial brightening, the first notable flash oc-curs at around 4s, followed by multiple bright flaring events.There are two distinct peaks in the light curve, at around5.5 s and at 5.9 s.We determined the total luminous energy of the fireballby integrating the light curve. We use only the visual magni-tudes (B, V and R). The corresponding total impact energy( E ) was calculated based on the total luminous energy ( E )using the empirical relation by Brown et al. (2002): E = . × ( E ) . (5)The impact energy then can be used, along with the pre-atmospheric speed to calculate the mass of the object. Theresults are shown in Table 4. The calculated average energyof the incoming bolide was about . ± . kt which yieldsan estimated mass of . ± . metric tonnes. − . − . − . − . − . − . − . A b s o l u t e M a g n i t ud e . . . . . . . Time [s] . . . . . L u m i n o s i t y [ k t / s ] UVBlueGreenRedIR
Figure 9.
The absolute magnitude and luminous power of thebolide from the different filters of the MSR as a function of time.The vertical dashed line indicates the peak brightness.
Table 4.
Results of light curve analysis. MSR E [ × J] 0.191[kt] 0.0455 E [ × J] 2.24[kt] 0.536Mass [metric tonne] 14.8
We determined the pre-atmospheric orbit from the velocityvalues measured at the earliest part of the bolide’s trajec-tory. First, we calculated the heliocentric position and veloc-ity vectors at the entry point. Then we carried out a back-ward time-integration using the REBOUND code (Rein &Liu 2012) to determine the position and velocity outside theHill sphere of the Earth. We tested several integration timesto ensure consistency in the retrieved orbit. There were smallvariations in the orbital parameters up to about 1.5 yearsbefore the event, after which the change was insignificant.Therefore, we integrated for 1.5 years, at which point we ob-tained the heliocentric orbital parameters. These values arelisted in Table 3.Figure 10 illustrates the orbit of the impactor whichwas inclined and slightly eccentric. The object reached per-ihelion about 64 days before the impact event. Figure 11
MNRAS000
MNRAS000 , 1–13 (2017)
Cs. Palotai et al. xy (a) yz (b)
Figure 10.
The diagram of the calculated orbit is shown in black,with the solid line above the ecliptic and dashed below in (a). (b)shows the view down the x-axis with the solid line before the y-zplane and dashed line, behind. The small and large black pointsshow the peri- and aphelion respectively. The orbits of Venus,Earth and Mars are shown along with their positions on June 2,2016. The vernal equinox is to the right in (a) and out of the pagein (b). shows the bolide’s orbital parameters plotted against the or-bital element distributions of known asteroids and comets.The object clearly did not belong to any known main-beltasteroidal families, but was part of the Apollo category ofnear-Earth objects.A search for the parent body was done using the Drum-mond criterion (Drummond 1981) and the modification byJopek (Jopek 1993) and data from the MPCORB database.This yielded no match for any known parent body. A significant number of movies became available on videosharing websites and news portals that were recorded bydash cams and security cameras. Most of these have theproblem of obstructed view throughout most or all of thebolide’s flight path and show only lens flare. In a few in-stances the entire event was clearly visible and for one ofthese videos we were able to identify the exact location whereit was recorded from. We have analyzed this video to giveus additional information on the bolide properties. Further-more, we have set up a fragmentation model to reveal piecesof the puzzle that was unavailable to us from the observa-tional data only.
Mark Olvaha recorded the event using his GoPro dash cam .Based on the land and road features seen in the video and inthe video description, we were able to pinpoint the locationof the car at the time of the recording as being on I-40 east h , accessedAugust 18, 2017 . . . . . . . . . . . . E cce n t r i c i t y . . . . . . Semi-major axis (AU) I n c li n a t i o n ( d e g ) Figure 11.
The asteroidal and cometary orbital element distribu-tion of known objects are shown in small black dots. The orbitalelement of the bolide is plotted with the large red circle. Thereare no distinguishable families that the bolide forms part of. of Kingsman, AZ at ( . N, . W). For our analy-sis we used a revised version of this video because in thisversion the authors applied image stabilization to the orig-inal source that made the bolide tracking easier. Since theposition where the dash cam movie was recorded was sig-nificantly further away compared to the SkySentinel cameralocations ( ∼ km), we were only able to get the light curveand the track of the object close to the end of the observa-tion (after t ∼ s). We calibrated the footage timing with theMSR dataset by matching the peaks of the sum pixel valuesfrom the video.Similar to the method described in section 2.4, the videowas calibrated for astrometry using local geographic featuresas shown in Fig 12. Features O, B, C, D, F, and G (six to-pographic points total) were used. Features A, E, H wereadditional markers that potentially could have been usedfor calibration but were not. A is a radio antenna that wasnot visible in the nighttime video; E is an interstate highwaysign whose azimuth changed significantly during the courseof the video; H is a mountain range without a well-definedsummit. The topographic features had well-defined silhou-ettes that facilitated determination of their azimuths and ap-parent elevations relative to the local horizon of the car. Us-ing Google Earth, it was straightforward to obtain azimuthfrom the known position of the car to each topographic fea-ture. To calculate the apparent elevation of each topographicfeature above the car’s horizon, a flat earth model was as-sumed. Each feature’s elevation above the car’s horizon wasobtained by first computing the difference in altitude (ex- - 0:00 to0:31, accessed August 18, 2017 MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Table 5.
List of geographical features used to calibrate the dash-cam footage. Feature AZ [deg] EL [deg]O 39.16 5.32B 48.81 5.29C 61.53 2.58D 80.76 3.09F 94.64 0.68G 115.55 1.95 pressed in terms of distance above mean sea level) betweenthe car and the feature’s summit. Each feature’s apparentelevation, el = arctan ( a / d ) was calculated, where a = altitudedifference and d = distance to the feature. Google Earth’stopographic data and measurement tools were essential tothis work.To generate an azimuth calibration model, the mea-sured x-coordinates of the six topographic features seen onthe stabilized video image were then regressed against theirazimuths. This resulted in an azimuth calibration modelwith R-sq (adj) = 99.9% and standard deviation = 0.98 de-gree. Similarly, the measured y-coordinates of the featureswere regressed against their apparent angular elevations.This yielded an elevation calibration model with R-sq (adj)= 98.0% and standard deviation = 0.26 degree. Minitab sta-tistical software was used to do the linear regressions.The calibration models had some limitations. The ap-parent elevation of the highest topographic feature was 5degrees, but the car’s video imagery first saw the bolidewhen it was about 15 degrees above the car’s local hori-zon. Extrapolation from 5 degrees to 15 degrees elevationnecessarily introduced significant uncertainty in the bolide’selevation. At an elevation of 15 degrees, the elevation un-certainty computed by Minitab statistical software is +/-2 degrees (95% confidence). The uncertainty of the bolide’sazimuth throughout its trajectory was also about +/- 2 de-grees (95% confidence). The AZ/EL of the features used aregiven in Table 5.By fixing the direction of travel (i.e. assuming that thebolide did not deviate significantly from this path), we wereable to convert the angular distance on the image to a lineardistance travelled by the bolide from the initial frame.This allowed us to constrain the distance travelled andthe speed of the bolide during and after the main fragmenta-tion, where the SkySentinel cameras were saturated or theirview were obstructed. The blue points in Figure 4 show theresults of the dash cam analysis. The error in the astrometryis on the order of about . ° which corresponds to a errorin the distance of about km. Therefore, we were unable toobtain a precise deceleration profile for the bolide. However,from the last few datapoints, it is clear that the bolide hasdecelerated significantly. A fit of the velocity for t (cid:38) s givesa velocity estimate of about - km/s. To obtain information on the deceleration of the bolide asobserved by the dash cam video, the position of the peakbrightness and the end of the luminous phase, we imple- mented the fragmentation model of Boroviˇcka et al. (2013)to calculate the post-breakup trajectory. We assume thatthe only processes that causes mass loss are discrete frag-mentation into one or more fragments and ablation.While erosion and the release of dust are important toaccurately calculate the mass loss, atmospheric trajectoryand the energy deposition, we opted to exclude these pro-cesses from the model because we could not properly validatetheir parameterization due to lack of precise velocity dataduring and post saturation of the bolide. To reduce the num-ber of free parameters, we also assume that only the mainbody breaks into smaller fragments and that ablation is theonly mass-loss process for the fragments.In this model, we input the times of discrete gross-fragmentation events (assumed to be local maxima in thelight curve), the mass loss during fragmentation and theestimated number of fragment produced. The initial pa-rameters for the model are the entry velocity, entry posi-tion, zenith distance and azimuth of the bolide radiant ob-tained from the SkySentinel triangulation. We approximatethe mass loss during each fragmentation event to be pro-portional to change in brightness of the respective peak inthe light curve, which allowed us to distribute the mass lossamong the flaring events. To be consistent with our energyanalysis, we use the luminous efficiency calculated using the(Brown et al. 2002) study, which in our case is . .Ceplecha & Revelle (2005) used results from observa-tions to constrain the values of K and σ . In our case, therewere no observations during the latter part of event, thus wecarried out sensitivity test for these unconstrained parame-ters. We varied the value of K from . to . (c.g.s) and σ from . - . s /km .The model is validated by matching the light curve aswell as the bolide positions from the SkySentinel and dash-cam video. The values of K and σ are kept constant through-out the trajectory, while the mass loss and number of frag-ments are varied to match the light curve. The end velocityof ∼ - km/s at t (cid:38) s, from the dashcam video, is used aconstraint to determine estimates of K and σ .In general, for K < . the bolide reaches an altitude of km with a velocity greater than km/s at t ∼ . s, whilefor K > . , there is insufficient energy to produce the brightpeak at . s. σ ∼ . s /km causes the bolide to brightenvery quickly at h (cid:38) km and mass loss to reduce the bolideinto sub-kilogram fragments well above km. For σ < . s /km the bolide continues to flare well after s. In eachof these cases, the mass loss at each fragmentation pointis kept constant. While these tests by no means produceunique solutions to the fragmentation model, they are usedto determine the type of object and minimize the number offree parameters to fit the observed light curve.To this end, we maintain the values of K = . cm /g / and σ = . s /km for the duration of the flight. The massloss during each fragmentation event is shown in Table 6,and resulting light curve is plotted in Fig 14. The resultingtotal luminous energy from the model is about 2% smallerthan the luminous energy from the MSR data.The bolide reached maximum brightness of − . ± . magnitudes at 5.5 s at a height of . ± . km. The endheight of the luminous phase is at a height of . ± . kmat . s at a planetographic latitude of . ± . ° N andlongitude . ± . ° W. MNRAS , 1–13 (2017) Cs. Palotai et al.
Figure 12.
Angular calibration of dashcam video using surface features.
Table 6.
Mass loss distribution during the fragmentation events.Time [s] Fragment mass [kg] Number of fragments4.00 7 204.50 27 4024 404.62 34 5025 454.91 18 4524 3527 305.35 14 5016 458 605.8 29 60
These parameters produced results that are in goodagreement with both the observed light curve and the tra-jectory. There are, however, a few instances where the modelpredictions differ from observational data:(i) The dashcam video shows that the bolide penetrateda few kilometers deeper than the fragmentation model pre-dicts, and the deceleration occurs at a much lower altitude.While it is possible to decrease the shape-density and ab-lation coefficients to match this, the resulting model bolidedoes not decelerate quickly enough to match the velocity at t ∼ . s. Therefore, it is more likely that inconsistency inthe penetration depth is an astrometric error due to the lowquality of the dashcam video.(ii) The peak at t ∼ . s from the fragmentation modelcontinues to flare for more than . s, which is contradictoryto the MSR data where the flare is very short-lived. Thisis most likely due to the fact that fragments were proba-bly smaller than the estimated by the fragmentation model.However, decreasing the modelled mass loss produces a muchdimmer peak, so either the bolide’s mass or velocity are inerror for the fragmentation event. Since there is no observa-tional data for the trajectory at this location, nor are there any observed fragments in the videos, it is not possible todetermine the exact set of parameters to match the lightcurve at this point.(iii) Initially the brightness of the model bolide is higherthan observed by the MSR data. We were able to reduce thisoffset by changing the ablation coefficient to . s /km at earlier parts of the trajectory, but it is unclear whythis would change throughout the trajectory to this degree.Since we have ignored the effects of dust release and erosionthroughout the flight, the contribution from these processeswould possibly explain the noted deviation.(iv) The peaks at t = . s and t = . s have a smooth flareas opposed to the sharp increase caused by the gross frag-mentation. These are likely eroding fragments as the massis slowly released from erosion. However, the absence of de-celeration information makes it impossible to constrain theexact mass loss from erosion and the size of grains (whichdetermines the intensity of the peaks), so we have insteadmodelled these peaks as the bolide exploding into numerous,small fragments. Based on analysis of weather Doppler radar images, re-searchers from the Arizona State University’s Center of Me-teorite Studies found fragments from the parent body onthe lands of the White Mountain Apache Tribe. 15 fusion-crusted stones were recovered during the 2 day search witha total of 79.46g of material (Garvie 2017). The locations ofthese Dishchii’bikoh meteorites are shown on Figure 3.
To calculate the size of the bolide from its mass, we requirean approximation for the bulk density value. We adaptedmultiple methods to determine this value indirectly fromobservations and we compare those results with properties ofthe actual “ground truth” value obtained from the recoveredfragments.
MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Table 7.
Result of size estimation from different methods. The errors in the sizes are carried through from the result for the total mass.Method Type/Material Bulk density [kg/m ] Diameter [m]PE criterion (Ceplecha & McCrosky 1976) IIIa 750 a . ± . Peak brightness (Brown et al. 2016) Type I 3000 a . ± . Type II 2000 a . ± . MSR visible bands analysis Fe-poor 2500 b . ± . Recovered Fragments (Garvie 2017) LL7 chondrite 3220 ± b . ± . a Ceplecha (1988) b Britt & Consolmagno (2003) V e l o c i t y [ k m / s ] Time [s] H e i g h t [ k m ] Figure 13.
The modeled velocity and height of the bolide areshown as a function of time. The solid black line is from thefragmentation model using the best fit values of K and σ , whilethe solutions of the triangulation are the red points. The solid redline in the velocity plot is the constant velocity approximationused in the triangulation from SkySentinel data. The blue box isthe approximate velocity constraint of − km/s for t between and . s from the dashcam video. Due to the uncertainty in thedashcam result, we were unable to determine the velocity withbetter precision. Figure 14. light curve from the model and MSR data.
Height [km] M a ss [ t o n ] Figure 15.
The mass of the main body as a function of heightfrom the fragmentation model. The vertical drops are from dis-crete fragmentation events, while the smooth mass loss is fromablation.MNRAS000
The mass of the main body as a function of heightfrom the fragmentation model. The vertical drops are from dis-crete fragmentation events, while the smooth mass loss is fromablation.MNRAS000 , 1–13 (2017) Cs. Palotai et al.
Ceplecha & McCrosky (1976) diagnosed the structureand bulk density values of impacting objects based on theterminal height of the fireballs they produced. Followingtheir method, the end height of km and entry velocity of . km/s points to a weak cometary material (Type IIIa).For a density of ρ ∼ kg/m , we get a size of . ± . m.However, such a weak material would probably not be ableto penetrate to a depth of ∼ km, so this is likely an inac-curate estimate of the bolide size.Brown et al. (2016) classified fireballs based on the anal-ysis of their peak brightness altitude values, using the samecategories. Using the appropriate parameters of the Ari-zona bolide, this object lies on the border of Type I/TypeII regimes. Using bulk density value of 2000 kg / m forType II objects (Ceplecha 1988), we obtained a diameterof . ± . m.For Type I objects, we take the density of about kg/m to be consistent with a weaker material as revealedby the fragmentation model. This gives a size of . ± . m. Looking at the SkySentinel MSR observations, the dataindicates strong emission in the blue filter and weak emissionin the green. This suggests that the object is Fe-poor, asa significant number of Fe emission lines are in the greenband. The high starting altitude of ∼ km points to alow material strength object of asteroidal origin as statedabove. For Fe-poor asteroids with low material strength, weadopted the bulk density value of / m was from Britt& Consolmagno (2003). Assuming spherical impactors, wecalculated an initial diameter of . ± . m for the bolide.Analysis of the recovered fragments categorize the ob-ject as being an LL7 chondrite with a shock stage S0 (Garvie2017). Britt & Consolmagno (2003) performed a bulk den-sity analysis for this material type and we adopted theirvalue of ± / m which results in an pre-impact di-ameter for the Arizona bolide of . ± . m. The values of K = . cm /g / and σ = . s /km , corresponding toa density of about - kg/m for Γ A (cid:48) ∼ . chosen inour fragmentation model are consistent with a weak objectof this type.The summary of results from all analyses describedabove is shown in Table 7. We have analyzed multi-station observations of a magnitude − . ± . superbolide that entered Earth’s atmosphere overArizona on June 2, 2016. The calculated deposited energyof the fireball from the light curve analysis is . ± . kt. This is about larger than the value of 0.48 kt re-ported by CNEOS based on data from US government sen-sors. Their reported peak location of . ° N and . ° Wlongitude (marked in Figure 3) is about 26 km off from ourcalculated position of the peak brightness. The site of recov-ered fragments are along the line of our calculated trajectory,thus we believe that for this particular event the CNEOS re-ported peak brightness location is inaccurate.From the deposited energy and the entry velocity of thebolide we estimate the mass of the object to be . ± . metric tonnes. We calculate the initial size of the object byassuming spherical shape. We use 3 different methods to ob- tain the bulk density of the meteor based on the observationsof the peak altitude, end hight, and spectral emission duringthe ablation. Boroviˇcka et al. (2017) noted that the methodof end height analysis (Ceplecha & McCrosky 1976) for su-perbolides may lead to misleading results, while the use ofthe height of peak brightness (Brown et al. 2016) to deter-mine material strength does not provides an adequate fit.Therefore, the true size of the object is most likely smallerthan what we obtained from these methods.Analysis of the MSR data points to a Fe-poor object. As-suming asteroidal origin, an Fe-poor composition, such asCI/CM chondrite, gives a slightly smaller diameter than thepeak/end height analysis. The reason we test different meth-ods is that in the absence of knowing the exact compositionwe would have to rely on these types of analyses. In thisparticular case, there were recovered fragments from thisbolide and using the actual LL7 chondrite composition weestimate the diameter of the object to be . ± . m. Ourspectral analysis is in agreement with this value while thepeak brightness method overestimate this diameter by about17%. In summary, the object is most likely to be ∼ m indiameter with a bulk density between − kg/m .The calculated orbit points to an origin interior to the mainbelt that is nearly co-orbital to the Earth. The object didnot form part of any known families, and a search for theparent body yielded no results. Dunn et al. (2013) state thata majority of Apollo-class asteroids are LL-chondrites; thisobject is therefore not unusual. ACKNOWLEDGEMENTS
The authors would like to thank J. Boroviˇcka for his thor-ough review and helpful comments that greatly improvedthe clarity of the manuscript. We also wish to acknowledgethose responsible for the operation and maintenance of theindividual Nodes used in the analysis of the Spalding AllskyCamera Network data.Nodes 1 and 2 - NMSU, Las Cruces, NM: Dr. Robert Wag-nerNode 6 - Flagstaff, AZ: Steven SchonerNode 14 - Los Alamos, NM: Dr. Matt HeavnerNode 37 - Parker, AZ: Jim Woodell (no longer operational)Node 79 - Turkey Springs Observatory, Payson, AZ: BruceRaschNode 7 - Lamy, NM: Dr. Thomas AshcraftNodes 5 and 8 - Albuquerque, NM: Dwayne FreeWe would also like to thank Dr. Jeremy Riousset for hishelpful suggestions.
REFERENCES
Astropy Collaboration et al., 2018, AJ, 156, 123Bannister S. M., Boucheron L. E., Voelz D. G., 2013, PASP, 125,1108Borovicka J., 1990, Bulletin of the Astronomical Institutes ofCzechoslovakia, 41, 391Borovicka J., Spurny P., Keclikova J., 1995, A&AS, 112, 173Boroviˇcka J., 1992, Publications of the Astronomical Institute ofthe Czechoslovak Academy of Sciences, 79Boroviˇcka J., et al., 2013, Meteoritics and Planetary Science, 48,1757 MNRAS , 1–13 (2017) nalysis of the June 2, 2016 bolide event over Arizona Boroviˇcka J., Spurn´y P., Grigore V. I., Svoreˇn J., 2017, Planetaryand Space Science, 143, 147Britt D. T., Consolmagno G. J., 2003, Meteoritics and PlanetaryScience, 38, 1161Brown P., Spalding R. E., ReVelle D. O., Tagliaferri E., WordenS. P., 2002, Nature, 420, 294Brown P., Wiegert P., Clark D., Tagliaferri E., 2016, Icarus, 266,96Ceplecha Z., 1987, Bulletin of the Astronomical Institutes ofCzechoslovakia, 38, 222Ceplecha Z., 1988, Bulletin of the Astronomical Institutes ofCzechoslovakia, 39, 221Ceplecha Z., McCrosky R. E., 1976, J. Geophys. Res., 81, 6257Ceplecha Z., Revelle D. O., 2005, Meteoritics and Planetary Sci-ence, 40, 35Ceplecha Z., Boroviˇcka J., Elford W. G., Revelle D. O., HawkesR. L., Porubˇcan V., ˇSimek M., 1998, Space Sci. Rev., 84, 327Drummond J. D., 1981, Icarus, 45, 545Dunn T. L., Burbine T. H., Bottke W. F., Clark J. P., 2013,Icarus, 222, 273Garvie L., 2017, Meteoritical Bulletin: Entry for Dishchii’bikoh,
Jopek T. J., 1993, Icarus, 106, 603Mirametrics Inc. 2017, Mira Pro x64,
Rein H., Liu S. F., 2012, A&A, 537, A128This paper has been typeset from a TEX/L A TEX file prepared bythe author.MNRAS000