A SPIRAL-BASED DOWNSCALING METHOD FOR GENERATING 30 M TIME SERIES IMAGE DATA

The spatial detail and updating frequency of land cover data are important factors influencing land surface dynamic monitoring applications in high spatial resolution scale. However, the fragmentized patches and seasonal variable of some land cover types (e. g. small crop field, wetland) make it labor-intensive and difficult in the generation of land cover data. Utilizing the high spatial resolution multi-temporal image data is a possible solution. Unfortunately, the spatial and temporal resolution of available remote sensing data like Landsat or MODIS datasets can hardly satisfy the minimum mapping unit and frequency of current land cover mapping / updating at the same time. The generation of high resolution time series may be a compromise to cover the shortage in land cover updating process. One of popular way is to downscale multi-temporal MODIS data with other high spatial resolution auxiliary data like Landsat. But the usual manner of downscaling pixel based on a window may lead to the underdetermined problem in heterogeneous area, result in the uncertainty of some high spatial resolution pixels. Therefore, the downscaled multi-temporal data can hardly reach high spatial resolution as Landsat data. A spiral based method was introduced to downscale low spatial and high temporal resolution image data to high spatial and high temporal resolution image data. By the way of searching the similar pixels around the adjacent region based on the spiral, the pixel set was made up in the adjacent region pixel by pixel. The underdetermined problem is prevented to a large extent from solving the linear system when adopting the pixel set constructed. With the help of ordinary least squares, the method inverted the endmember values of linear system. The high spatial resolution image was reconstructed on the basis of high spatial resolution class map and the endmember values band by band. Then, the high spatial resolution time series was formed with these high spatial resolution images image by image. Simulated experiment and remote sensing image downscaling experiment were conducted. In simulated experiment, the 30 meters class map dataset Globeland30 was adopted to investigate the effect on avoid the underdetermined problem in downscaling procedure and a comparison between spiral and window was conducted. Further, the MODIS NDVI and Landsat image data was adopted to generate the 30m time series NDVI in remote sensing image downscaling experiment. Simulated experiment results showed that the proposed method had a robust performance in downscaling pixel in heterogeneous region and indicated that it was superior to the traditional window-based methods. The high resolution time series generated may be a benefit to the mapping and updating of land cover data. * Corresponding author


INTRODUCTION
Remote sensing image with high spatial and high temporal resolution acts as a critical role in environmental change monitoring, land cover data updating and land resource management (Chen et al. 2015).However, due to many factors such as satellite orbit, sensors hardware and so on, it is difficult to provide high spatial and temporal remote sensing image in large area.MODIS sensor in Terra and Aqua repeats global observation four times within one day, but its highest image spatial resolution is only 250 meters.One popular way is utilizing the texture information in sparse high spatial resolution image and phonological information in dense high temporal resolution images to generate corresponding high spatial resolution image series (Gao et al. 2006).With some assumptions, the linear mixing spectrum decomposition based approaches downscale coarse pixel to fine endmember by taking into account adjacent coarse pixels and solving its linear system (Wu et al. 2015).Unfortunately, the endmember types are usually different in coarse pixel adjacent region in heterogeneous area, leading to the unmixing linear system have no unique solution.As a result, some endmembers have unreasonable value out of bounds.
Usually, the adjacent region is defined by a regular window with certain size, such as a rectangle region of 3×3, 5 ×5 pixel size center on a coarse pixel.This way is mostly effective in land surface, but its regular adjacent region can't portray the land cover exactly.As a result, the linear system may lose some endmember types ensuring correct calculation.In other words, it brings the uncertainty in the downscaling of coarse pixels (He and Gaston 2007).Researchers weight the pixels among the regular region demarcated by window (Busetto et al. 2008).These improvements raise the overall accuracy to some extent, but can't avoid underdetermined problem.Others set different window size based on land surface heterogeneity (Zhu et al. 2010).However, this approach is still not the fundamental solution.One of nature factors is regular rectangular windows are hard to accurately approximate the spatial distribution patterns of the objects on the land surface.For example, crop, water and artificial land cover type distribute irregularly along the river region.It is difficult to conform to the trend of the river for regular window.Accordingly, the linear system can't get enough proper endmember types to be solved successfully.
The spiral is a kind of geometric curve with irregular spatial form.Compared with regular window, spiral is closer to the irregular land surface around the pixel adjacent region.Along with the spiral path, it is easy to select pixels with similar endmember types to construct pixel set around adjacent region.When the pixel set contains enough proper pixels, the underdetermined problem could be solved fundamentally.

Spiral-based solution
The Arithmetic spiral is a curve that winds around a fixed central point with a constant speed along a line while rotating with constant angular velocity.A spiral center on the origin can be described by the equation (1) in polar coordinate, which refers to polar radius and θ refers to polar angle，as shown in figure 1.Compared with regular window, spiral is more flexible to describe the irregular land surface around the pixel adjacent region.Along the spiral path, it is easily to select pixels with similar endmember types to construct pixel set around adjacent region.When the pixel set contains enough proper pixels, the underdetermined problem could be solved fundamentally.

Figure 1. A spiral in polar form
In figure 2, the class map contains crop, forest, grass, shrub, wetland, water, artificial and bareland, total land type numbers is 8. Four spirals with different length are adopted to construct adjacent pixel set for coarse pixel A, B, C and D. Letter A in figure 2 stands for case that pixel A contains forest, grass and crop types , the spiral select other four pixels with similar endmember types; letter B in figure 2 stands for case that pixel B contains crop, water and wetland, the spiral select other two pixel with similar endmember types; letter C in figure 2 stands for case that pixel C contains crop, bareland, artificial, forest and grass, the spiral select other four pixel with similar endmember types; letter D in figure 2 stands for case that pixel D only contains crop, the spiral select pixel D as adjacent pixel set.The type and distribution of land objects is various in land surface.The adjacent objects relate each other strongly in landscape homogeneous region and weakly in landscape heterogeneous region.For coarse pixels, it means that the inner endmember types are similar among adjacent region; otherwise not.It is desirable to take a elastic manner relative to land surface facing for various coarse pixels.Compared with window, the spiral is more adapted in adjacent region.

Coarse pixel downscaling
For the remote sensing image, it is a mixed pixel decomposition process from low spatial resolution scale to high spatial resolution.In order to decompose the coarse scale pixels into the fine scale, the mapping relation between the coarse and fine scale pixels is firstly fitted.The linear spectral mixture model (LSM), because of its physical meaning of parameters, has been widely used.
In the linear spectral mixture model, each pixel value is expressed as a linear combination of the components and the proportion of the components in the pixel.In general, the formula of the linear spectral mixture model can be formulated as equation ( 2), where R is the low spatial resolution scale remote sensing image pixel value, f is component proportion in the low spatial resolution range, n pixel components in total, error is nonlinear components other than linear mapping relation in the model.With this model, each pixel in the spiral based pixel set was translated into one equation.All equations made up the linear system, which the solutions stand for endmember values.
After the classification of Landsat image, the class map with 30m spatial resolution was got to drive the abundance of endmembers within each coarse pixel.Then, these abundance were used to invert the endmember values.The class map provided the high spatial information and the inverted values provided the high temporal information.Based on the class map, these inverted values were constructed to form the high spatial resolution image band.The high spatial resolution time series was made up band by band and image by image.

Study area and data
In order to verify the performance of the proposed method, an experiment based on GlobeLand30 dataset was conducted.The study area covers about 1300 km2 land surface and is located in the east longitude between 115.6 °-117.3°and north latitude between 39.4°-40.5°ofBeijing area in China.Forest and grass distribute widely in northwestern mountain area.The fragmented artificial landscape are mainly located in southeastern plain area.Crop land accounted for the majority and some rivers and lakes are among them.Overall, land surface types are crop, forest, grass, shrub, wetland, water, artificial surface and bare land.This study area is representative of the landscape with heterogeneity and homogeneity.

Results
The simulate experiment count the total underdetermined problem pixels within study area by proposed method and present with an unit of 1/10000.By constructing adjacent pixel set on low spatial resolution, each pixel to be unmixed builds a linear equation system solved later.It is marked as underdetermined problem pixel when one pixel corresponding system cannot be solved.A comparison between spiral-based and window-based method results is made respectively and discussed the applicability of them.For linear mixing spectrum decomposition, the popular approach is unmixing the low spatial resolution (250 m) MODIS NDVI pixel to high spatial resolution (30 m) Landsat NDVI liked pixel with according land cover map.For our land cover map, this means that one low spatial resolution pixel occupies about 8 × 8 high spatial resolution pixels.
In figure 4, the underdetermined problem pixels are marked as black.From figure 4 a and e ,we can find that much grass patch locate in forest edge.That results in complicated land cover type constitute within a small region, especially in a adjacent region.The window-based method can hardly handle this situation, the black pixel mainly scattered in these areas.As a contrast, the spiral-based method deals with this situation better.From figure 4 b and f, it can be found that some black pixels located in artificial land edge in window-based method and almost none in spiral-based method.From figure 4 c and g, we can find that a river located in crop and artificial land.Some small clustered artificial land patch and dotted water patch intertwined along the river region.As a result, it increases the difficulty of decomposition.The spiral-based method adapts this situation well.In figure 4 d and h, fragmentized land cover patch scattered along lake.The complicated combination of grass ,forest, crop and water bring great difficulties for constructing adjacent pixel set.In this case, the spiral-based method works better as well.Among the forest, water, artificial land and crop land areas, the black color pixels on spiral-based method results is about 0.9‱ .In the same region, the black pixels on window-based method results is about 42.7‱ .
In the city area, the real NDVI data of downscaling results as shown in figure 5.By comparing figure 5 about two images, we can see the whole downscaling effect and is close to the real NDVI data, better integration.

CONCLUSION
This research revealed the imperfection of traditional window based methods in linear mixing spectrum decomposition and proposed a novel downscaling method based on spiral.The spiral based method constructs the adjacent pixel set along with a rasterization spiral route and selects the similar pixels.The experiment demonstrated that this new method could achieve a robust capability for homogeneous or heterogeneous region and improved the underdetermined problem pixels in the linear system to a great extent.The downscaled NDVI results had a low error compared with the true NDVI data.
However, the spiral based method still has constraints.For example, when the mixed pixel contained two unique endmember types among the whole research region, the method could not construct the adjacent pixel set to avoid the underdetermined problem all the same.The proposed spiral based method will reduce the hard unmixed pixels caused by the underdetermined problem in linear mixing spectrum decomposition greatly and improve the accuracy in other unmixing / downscaling researches by integrating with this spiral based method.

Figure 2 .
Figure 2. Examples of spiral-based method constructs adjacent pixel setThe type and distribution of land objects is various in land surface.The adjacent objects relate each other strongly in landscape homogeneous region and weakly in landscape heterogeneous region.For coarse pixels, it means that the inner endmember types are similar among adjacent region; otherwise not.It is desirable to take a elastic manner relative to land surface facing for various coarse pixels.Compared with window, the spiral is more adapted in adjacent region.2.2 Coarse pixel downscaling

Figure 4 .
Figure 4.The underdetermined pixels among class map(fig.a, b, c, d were window based results, fig.e, f, g, h were window based results)