Received 9 February 2016; accepted 14 May 2016; published 17 May 2016
Sunlight passing through the atmosphere is scattered and absorbed by gases and aerosols to produce the characteristically nebulous appearance. The amount of scattering and absorption that light undergoes depends upon the abundance of aerosols as well as the wavelength of the radiation. Light at certain frequencies interacts with the gaseous molecules and aerosols more effectively than at others; consequently, the depth to which light effectively penetrates and escapes the atmosphere is a function of wavelength. By using this wavelength dependence, it is possible to image clouds at different altitudes. Clouds, depending from both their cold/hot nature and their growing mechanism, can assume a variety of shapes and sizes. In this work, we propose a method suitable to estimate clouds quite elliptically shaped. Though this technique could be applied to any planetary cloud system, it has been setup aiming to infer the evolution of a particular vortex observed on Saturn in occasion of a giant storm sprouted at the end of 2010. Ground-based telescopes detected the first signs of the storm on 5 December 2010 at 21:05:22 Universal Time (UT), but it was possible also to observe and probe it by the Cassini mission. Cassini Composite Infrared Spectrometer (CIRS)  , Radio Plasma Wave Science (RPWS)  , Imaging Science Subsystem (ISS)  and Visual and Infrared Mapping Spectrometer (VIMS)  sounded different pressure levels and exploited a variety of features associated to the storm development. Many works have already been carried out on this topic by means of the CIRS thermal infrared observations   , the RPWS radio pulses  and the camera images analysis  in the 2011 January to October lapse time.
The storm structure, in upper atmosphere, can be usefully analyzed by ISS and VIMS. The instruments cover the UV (only ISS), the visible (VIS) and the Near Infra Red (NIR) spectral ranges. They observed the giant storm on Saturn during its whole life cycle, with some regularity. The first detections of the outbursts’ tail started on February 2011, but the feature object of this report has been observed, in various occasions and from different distances, until December 2012. Here, we will describe the position and morphology history of that feature from its first occurrence up to January 2012. We assume that the vortex monitored by VIMS is the same described in  , for its position, similar shape and remarkable dimensions, thus we analyze the series of measurements reported in Table 1, starting with the ISS observation on 12 January 2011 and ending with the VIMS observation on 4 January 2012. In Section 2, we will describe the data selected for the analysis and the processing algorithms that produce the required outcomes. In Section 3, we will discuss our results and summarize our work.
2. Observations and Data Handling
In this paper, the observations by the VIMS visual channel and by ISS-WAC (Wide Angle Camera), listed in Table 1, are analyzed. VIMS is a hyperspectral imaging spectrometer; its visual channel (hereafter VIMS-V) has 96 bands ranging from 350 to 1046 nm with a 7 nm nominal bandwidth. ISS-WAC is equipped with two filter wheels with a wide variety of filters. Both VIMS-V and ISS-WAC raw data are publicly available from the NASA Planetary Data System (PDS) Imaging Node. Our search on VIMS-V datasets has been carried out principally on the archive of IAPS-INAF, where updated geometrically and radiometrically calibrated data are available for the Italian scientific team. ISS observations were selected with the purposes of filling the large gaps in time in which VIMS-V was used to acquire data for other scientific objectives, or when Cassini’s orbit was not adequate to observe the relevant latitude. In Table 1 the list of the VIMS-V and ISS-WAC observations is reported together with the respective pixel sizes, obtained through the conversion of the images from a geographical to a Universal Transverse Mercator (UTM) projection.
The projection calculations have been done exploiting some tools and facilities from ENVI-IDL (www.exelisvis.com) and ESRI (www.esri.com) software. UTM projection provides a north-south oriented swath of little distortion. Each of these swaths is called a UTM zone and is 6˚ wide in longitude. The first zone
Table 1. VIMS-V and ISS-WAC observations used for the vortex study.
conventionally begins at the International Date Line (180˚, using the geographic coordinate system). Saturn_2000_IAU_IAG ellipsoid and D_Saturn_2000 datum are the ESRI’s approximation of the Saturn’s irregular shape and the mathematical model used to describe the location of unknown points on Saturn’s surface, respectively. This geographic information, along with the pixel geo-referencing by C-SPICE kernels  , permitted us to evaluate the spatial frames’ dimension and orientation. The raw VIMS-V data are radiometrically calibrated by using the official calibration pipeline   . The VIMS and ISS data have been geo-referenced using ad hoc algorithms based on NAIF-SPICE tools. No radiometric calibration has been applied to ISS data since in this paper we only focused on the morphology (shape and size) evolution of the vortex and on its displacement in the examined time period. The ISS images inspection aimed principally to maximize the probability of a univocal recognition for the structure under analysis, adopting the same criterion used in  . We included two ISS-WAC measurements in our survey: on 12 January, a time period without VIMS-V observations, and on 24 August, when both instruments, in about half an hour and with the same observation geometry approximately, observed the same region. The 24 August observations have been used to compare the intensity scales of the two instruments and to check their correspondence. The two images have been resized on a geographic base, re-sam-pled to a common dimension of pixel number and normalized on a monochromatic color range scale (8 bit, 0 - 255 levels). In the resulting images (Figure 1) the contour lines of the prevalent structures visible at 750 nm of VIMS (left panel) and of ISS (CB2 filter, right panel) are equivalent in spite of the different spatial resolutions. We are confident that, if in this complex dynamic case similar values have been obtained, an analogous behavior can be assumed also for other wavelengths and simpler configurations.
3.1. Size, Shape and Dimension Evaluation
The evaluation of the limits of the feature under study presents the problem of identifying its boundary without ambiguities. Not all the spectral images at different wavelengths show the same appearance of the vortex. As can be seen in Figure 2, where three spectral frames of the January 2012 VIMS-V observation are shown, there is a great difference between the vortex shape as revealed at the ultraviolet wavelength (left panel, 373 nm), in correspondence of a methane window (central panel, 752 nm) and in the deepest methane band of the visible range (right panel, 890 nm). Each one of these wavelengths sounds the atmosphere at a different depth. The UV wavelength (left panel, ~400 nm) senses relatively deeply, but the color scale (corresponding to spectral radiance on the sensor) is influenced mostly from particles’ absorption and Rayleigh scattering by gases. The 752 nm (central panel, high continuum reflectivity) sounds the deepest level and is the least affected by gas scattering. In the methane absorption band (right panel, ~890 nm) the reflected light comes principally from the highest particles, because methane absorbs both the transmitted and the diffused radiation deepen in atmosphere. Thus the bright oval feature in Figure 2(c) corresponds to a high altitude light reflection, well-contrasted against a darker, deeper background caused by methane absorption. This study is carried out on the third case because it offers a simpler structure showing the top vortex as a bright oval (BO in the following). However a quantitative criterion to choose where inserting the boundary line between the BO shape and the background had to be fixed also in this case.
The VIMS-V and the ISS-WAC data have been normalized in order to make both datasets homogeneous by:
where S is the signal (spectral radiance for VIMS-V and Digital Number (DN) for ISS-WAC), Smin and Smax the minimum and maximum values for each dataset, and SN is the normalized value, respectively.
In Figure 3, the spectral frames at 890 nm for the observations listed in Table 1 before (top) and after (bottom) the re-sampling and normalization are reported. The nearly elliptical shape of BO at 890 nm has then been fitted with a two-dimensional, elliptical Gauss function in order to objectively determine its boundaries, sizes and shapes. The fitting equation F (x, y) = A0 + A1 × exp (-U/2) has been applied to the observations reported in Table 1, but for COISS_2069-W1692857065_1, referring to 752 nm. In the equation U is the elliptical function U= (x’/a)2 + (y’/b)2, and A0 and A1 are the intercept and the scale factor of the elliptical Gauss function. In the U equation a and b are the FWHM of the Gaussian in the X and Y direction (sample and line in the pixel coordinate system), respectively. A code from the IDL library has been used for this calculation. This procedure requires to operate on a rectilinearly gridded data (Figure 4). The fitting routine is based on an adaptation from
Figure 1. In the two panels the vortex observation on 24 August 2011 is reported in a monochromatic color range scale (0 - 255 levels). Both images have been resized through the associated geographical coordinate system (GCS). Left panel: the VIMS-V observation acquired on 04:55 UTC mean time. The image corresponds to the band 56 centered at 753.96 nm. The spatial frame has been re-sampled to the ISS-WAC pixel number dimensions. Right panel: the ISS-WAC image acquired on the image-mid-time 05:15 UTC through the CB2 filter (752 nm).
Figure 2. The vortex observed by VIMS-V on 4 January 2012 at three different wavelengths. Light imbalance has not been corrected. Sun illumination and viewing angles grow from right side to left side of the image: (a) The structure at 373 nm appears as an elliptical shape encircled by a darker ring; another smaller oval is visible on the bottom-right corner of the image. (b) At 752 nm the vortex shows many well-contrasted dynamic details; the previous smaller oval corresponds to the brightest feature at this wavelength. (c) The deepest methane band at 890 nm enhances the high feature, located on the top of the vortex; the smaller oval appears darker (deeper) than its background.
Figure 3. Top: spectral frames at 890 nm by ISS (the first at left) and VIMS-V observations. At this wavelength the CH4 layer absorbs the radiation going deep in the atmosphere and only the light scattered by the stratospheric hazes is visible. The different spatial resolutions due to observation distances and device characteristics make difficult to fix a size for comparison. Bottom: normalized and re-gridded values with contour plots.
Figure 4. Top: in yellow the images of Figure 3 are visualized as surfaces in the X (sample), Y (line), Z (normalized units) coordinate system together with the overplotted best fitting elliptical Gauss function (green) and the contour lines of the surface projected on Z = 0 plane. Bottom: the images of Figure 3 after the pre-processing and with the fitting ellipses overlaid.
Marquardt  , elaborated by Bevington and Robinson  and returns the χ2 value, as goodness of fit parameter, a seven element vector containing, among the others, the widths of the Gaussians in the X and Y direction, the ellipse center location in pixel coordinate system, and, if present, the tilt angle of the ellipse in the clockwise direction from the X axis in radians.
In Figure 4 (top), the best fitting surface (green) over-plotted to the surface representation of BO (yellow), imaged in Figure 3, is shown. Contour of the normalized values are plotted on the Z = 0 plane. At the bottom of Figure 4, the fitting ellipses are overlaid to the rectilinearly gridded images of Figure 3.
The U ellipse parameters actually refer to X (sample) and Y (line) couples in pixel coordinate system. Thus couples of longitude and latitude values can be associated to each geo-located pixel and corresponding am, bm pathlengths can be calculated for the fitting ellipse in metric units (see the eccentricity equation ahead). This operation implies a deformation of the original shape in pixel coordinates, connected to the warping that image undergoes once projected. The fitting ellipse axes have been converted in their corresponding pathlengths through another code of the IDL library that returns the rhumb line-the line crossing all meridians of longitude at the same angle-that connects two points on a sphere, known the latitude and longitude of the two points and the radius of the sphere. As Saturn is not a spherical body, the local radius of curvature of the meridian passing for the vortex center had to be calculated in place of the sphere radius. The BO shape has been estimated through the eccentricity value in pathlengths metric units, or:
where am and bm are the major and minor half pathlengths in metric units. The size corresponds to the U area values. The shape tilt on the equator line is calculated as the regression coefficient of the ellipse major axis, or on the three couples of longitude (x), latitude (y) values of the U center and the two vertices of the largest pathlength. In Figure 5, the U ellipses, imaged in Figure 4 (bottom) in image plane coordinate system, are visualized on the map of the observations in planetocentric coordinate system. Size, shape, areas and tilt are given in Table 2.
3.2. Error Analysis
Errors on the BO size, shape and dimension are derived from the kilometric pixel size reported in Table 1. The kilometric pixel error is calculated as the diagonal half value of the square region around the pixel center. Errors connected to the fit parameters are negligible compared to the size error.
The eccentricity (shape) error is obtained from the propagation of the pathlengths uncertainties.
Analogous computation has been carried out to evaluate the error on the area.
3.3. Drift Evaluation
Latitude and longitude drift has been estimated following the BO center coordinate displacement vs. time. Planetocentric longitudes have been offset to start from the 0˚ meridian as initial value. Really five points have to be considered with some caution to infer a trend. However, considering the quite regular behavior and the typical trends from previous works on similar topics   , some conclusions can be draw. Then we have interpolated
Figure 5. The same of Figure 4 but in planetocentric coordinate system.
Table 2. χ2 value, center location, axes size, shape, area and tilt for the BO observations reported in Table 1.
Figure 6. Displacement and propagation rate of BO center are illustrated for the meridional (upper panel) and zonal (lower panel) components. Center coordinates and the best fit curve (red) values refer to right axis. First derivative curve (black dashed lines) refers to left axis. In abscissa the days since the first observation on 12-01-2011 are reported.
the five points of the ellipse centers in latitude and longitude with polynomials fit, where the x axis is a Julian day’s grid starting from the first observation of Table 1. The first derivatives of these curves have been also calculated on each grid point to obtain a propagation rate trend and finally the average on the whole time interval has been assumed as drift speed value.
BO results to have a negligible mean propagation rate of 7.16 × 10−3 ± 6.66 × 10−3 deg/day toward the North. The mean propagation rate in westward direction is instead much faster being 0.75 ± 0.25 deg/day, corresponding to 8.9 ± 3.0 m・s−1 in System III reference frame, in agreement with the results in  . It is worth to observe
Figure 7. The area values of the fitting ellipses indicate a shrinking with the time of the BO in agreement with Fletcher et al. (2012) and Sayanagi et al. (2013) studies.
that the only ISS-WAC observation used in this report is not decisive in the drift speed determination, however confirming the tendency of the BO to propagate in clockwise direction. In Figure 6, the BO center coordinates (red points), the best fit of the displacement (red line) and the first derivative of the best fit displacement (black dashed line), are reported both for the latitude (top frame) and for the longitude (bottom frame).The average propagation rates for the meridional (top) and zonal (bottom) components are depicted as black solid lines. Watching at the first derivative behavior, the BO meridional displacement would seem undergoes an acceleration with the time, but the low values associated (left y axis) induce to consider it stationary in latitude (top frame). Zonal acceleration of BO (bottom frame) shows instead a fast raising in the first 100 days to reach a maximum around the day 140, followed by a slow decrease. This is in line with a dynamic connected to the outburst interpretation, where the strong initial energy of the diffusive wake is damping with the time.
Summarizing the results reported in Table 2, the BO area seems to undergo a shrinking up to January 2012 (Figure 7). The shape evolution does not show a similar linear behavior, oscillating apparently in a random way during the first half of the year, as well as the tilt, and approaching a circular form at the end of the period. On the other hand the fit has been done on the shape of the image plane, thus a sensible uncertainty on the effective contour have to be considered as a bias.
In conclusion, the mathematical procedure setup to analyze the life cycle of the top cloud over a vortex in the Saturn’s giant storm on 2010 has been successfully checked. Dimensions are very different from those reported in  , however in this work different wavelength frames and pressure levels of the cloud features are considered. On the other hand, size shrinking, direction and speed drift are confirmed. Future work aims to apply this technique to other spectral frames with the purpose to produce a 3D model of the vertical structure, which would help to constrain the related dynamic origins.
Special thanks are given to the Cassini/VIMS and/ISS teams for planning the observations. The ISAC-CNR and IAPS-INAF teams have been funded by the Cassini project of the Italian Space Agency under the contract ASI I/015/09.
 Flasar, F.M., et al. (2005) Titan’s Atmospheric Temperatures, Winds, and Composition. Science, 308, 975-978.
Acton, C.H. (1996) Ancillary Data Services of NASA’s Navigation and Ancillary Information Facility. Planetary and Space Science, 44, 65-70.
 Porco, C.C., et al. (2004) Cassini Imaging Science: Instrument Characteristics and Anticipated Scientific Investigations at Saturn. Space Science Reviews, 115, 363-497.
 Fischer, G., Kurth, W.S., Gurnett, D.A. Zarka, P., Dyudina, U.A., Ingersoll, A.P., Ewald, S.P., Porco, C.C., Wesley, A., Go, C. and Delcroix, M. (2011) A Giant Thunderstorm on Saturn. Nature, 475, 75-77.
 Sayanagi, K.M., Dyudina, U.A., Ewald, S.P., Fischer, G., Ingersoll, A.P., Kurth, W.S., Muro, G.D., Porco, C.C. and West, R.A. (2013) Dynamics of Saturn’s Great Storm of 2010-2011 from Cassini ISS and RPWS. Icarus, 223, 460-478.
 Filacchione, G., et al. (2007) Saturn’s Icy Satellites Investigated by Cassini-VIMS. I Full-Disk Properties: 350-5100 nm Reflectance Spectra and Phase Curves. Icarus, 186, 259-290.
 del Río-Gaztelurrutia, T., Legarreta, J., Hueso, R., Pérez-Hoyos, S. and Sánchez-Lavega, A. (2010) A Long-Lived Cyclone in Saturn’s Atmosphere: Observations and Models. Icarus, 209, 665-681.
 Vasavada, A.R., Hoerst, S.M., Kennedy, M.R., Ingersoll, A.P., Porco, C.C., Del Genio, A.D. and West R.A. (2006) Cassini Imaging of Saturn: Southern Hemisphere Winds and Vortices. Journal of Geophysical Research, 111, E05004.