ASSESSING THE CONTRIBUTION OF SPECTRAL AND TEMPORAL FEATURES FOR ANNUAL LAND COVER AND CROP TYPE MAPPING

Freely available satellite image time-series are currently the most exploited data towards land cover mapping. In this work we assess the contribution of spectral and temporal features for the detailed, i.e., with more than thirty classes, land cover and crop type mapping based on annual Sentinel-2 data. As a baseline we employed a datacube consisting of spectral features, i.e., spectral bands and indices from one tile of Sentinel-2A data for the year 2016. Then we formed two different datacubes of reduced dimensions, containing either spectrotemporal or temporal features and performed the same experiments in order to assess their contribution. For the second dataset only spectral features that fulfilled certain temporal conditions were retained, reducing by 40% the initial datacube dimensionality. The third dataset was formed only of temporal features resulting to a reduction of 50%. A random forest classifier was employed for the classification procedure and standard accuracy metrics for the validation. All experiments resulted into high overall accuracy rates of over 90% while rates for average F-score metric exceeded 78% in all cases. The quantitative and qualitative validation indicated that the baseline dataset modestly outperformed the other two of spectrotemporal and temporal features. Insights regarding the influence of spectral differentiation among classes and the impact of their sample size, on the per-class performance are further discussed. The importance of spatial independency for training and testing sets was also demonstrated highlighting the need of following best practises during validation in order to deliver a realistic estimation of the produced map accuracy.


INTRODUCTION
Accurate and regularly-updated mapping along with change analysis arises as essential for several scientific communities, but also for public and regional authorities in terms of supporting decision making, planning, sustainable development and natural resources management. At the same time open data policies both in the USA and EU, are delivering an unprecedented volume of satellite data with increasing levels of resolution. Currently, the availability of Landsat-8 and Sentinel-2 data significantly expands the capabilities of timely, accurate and detailed land cover mapping from time series of cost-free satellite data. The use of multi-temporal data for land cover mapping tasks at annual basis has become currently the standard method, since in this way the valuable information of the phenological variations of different land-cover types can be exploited (Cihlar, 2000;Inglada et al., 2017;Jia et al., 2014;Karakizi et al., 2018b;Xie et al., 2019) Mapping crop areas and classifying different crop types, arises as an even more complex problem requiring the use of data with higher spatial, spectral and temporal resolution. In recent studies the demanding task of crop-type mapping has been addressed with the use of satellite imagery combined with machine learning frameworks (Defourny et al., 2019;Karakizi et al., 2018a;Lira Melo de Oliveira Santos et al., 2019). The use of machine learning techniques like Support Vector Machines (SVM) and Random Forests (RF) have gained rapid recognition for land cover and crop type classification studies (Defourny et al., 2019;Inglada et al., 2017Inglada et al., , 2015Karakizi et al., 2018bKarakizi et al., , 2018aXie et al., 2019;Zhai et al., 2018). A different approach, deep learning and Neural Networks (NN), is currently one of the fastest-growing trends in remote sensing. Compared to more shallow architectures (like SVM and RF) for the task of land cover mapping, the NN based analysis has given results of similar accuracy. However, it is accompanied by a significantly increased computational cost related with the complexity of the training procedure combined with the high dimensionality of image time series (Karakizi et al., 2018b;Khatami et al., 2016;Stoian et al., 2019).
Apart from the optimal selection of a classifier that balances highly accurate results and lower computational needs, the choice of classification features also plays an important role in the efficiency of the classification framework. Spectral features, namely multispectral satellite bands and derived spectral indices, have been established as the main set of input features for land cover and crop type classification in the recent literature (Defourny et al., 2019;Inglada et al., 2017Inglada et al., , 2015Karakizi et al., 2018bKarakizi et al., , 2018a. High dimensionality issues and redundancy problems of unchanged regions related with multispectral time series are often tackled by deriving temporal metrics from the multi-temporal data (Zhai et al., 2018). Multitemporal metrics are statistical derivations of the time series per pixel, e.g., maximum, minimum, median and other percentile values of the spectral bands.
Temporal features are frequently used and recommended for crop-areas mapping (Defourny et al., 2019;Song et al., 2017;Valero et al., 2016;Waldner et al., 2017) since those type of metrics capture vegetation phenology but are generally insensitive to the timing of phenological differences across large areas (Egorov et al., 2018). Temporal metrics can also be employed to produce spectrotemporal features. This is achieved by defining key dates of remote sensing stages for the dynamic classes (vegetation, crops etc.) and then deriving the reflectance values of the spectral bands on those specific dates. As a result, the dataset to be classified is independent of the calendar sequence, i.e., the sequence of acquisition dates, for the intraannual changes and thus grants better handling of the dynamic The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII- B3-2020, 2020XXIV ISPRS Congress (2020 This contribution has been peer-reviewed.
cropland diversity and the agro-climatic gradient across the landscape (Matton et al., 2015).
Towards this direction, in this paper we assess the contribution of spectral and temporal features for the challenging task of detailed land cover and crop type mapping based on multitemporal Sentinel-2 data from the year 2016. Aiming to reduce the computational costs of a hardware-hungry baseline dataset consisting of spectral features, we formed two additional datasets with reduced dimensions. The reduced datasets were formed employing spectrοtemporal and temporal features derived from the initial baseline dataset. The same experimental set-up based on a RF classifier was applied to all three datasets, to ensure an objective comparison. A comprehensive quantitative and qualitative evaluation was performed for all experiments, while additional aspects like the spatial independency between training and testing sets were also discussed in order to deliver a reliable estimation of the produced map accuracy.

Study Area
The study area is located in central continental Greece and corresponds to the Sentinel-2 tile 34SEJ ( Figure 1). It covers an extent of about 10,000 km 2 and includes parts of four administrative regions of the country. We selected this region as the study area since it presents highly heterogeneous landscapes and varying terrain relief including plains but also several mountain masses. The Pindus mountain range crosses the study area from north-west to south-east, presenting covers of natural vegetation such as broadleaved and coniferous forests, natural grasslands, sclerophyllous vegetation and barelands. On the east of Pindus lies a big part of the Thessalian Plain, one of the most important agricultural regions of Greece. Agricultural land consists mainly of cotton, maize, cereals, clover and other grass fodders. Urban areas like the cities of Trikala and Karditsa and smaller towns and villages are scattered across the plain. The study area also includes several water bodies, e.g, man-made lakes Plastira, Pournariou, Kremaston, while on the south-west border of the study area, lies the salty-water lagoon of Messolonghi. Figure 1. The study area (in red box) corresponding to the Sentinel-2 tile 34SEJ (in a natural color composite of March 2016) and its location on the Greek territory

Sentinel-2 Data Selection and Preprocessing
Sentinel-2A L1C data of less than 10% cloud coverage (10 dates) for year 2016 was downloaded from the ESA Sci-Hub for the tile 34SEJ. An atmospheric correction procedure followed using the Sen2Cor processor, in order to produce surface reflectance products of level 2C. Afterwards, BRDF correction was carried out, where pixel reflectance values were adjusted as if the satellite's location was at the nadir concurrently with image acquisition. Then, cloud and shadow screening was conducted with the F-mask algorithm and interpolated values using the previous and following cloud/shadow-free dates, were produced for cloudy and invalid pixels. Furthermore, geolocation shift errors (1-2 pixels) that may occur between different dates/images of the same tile were also resolved. As a final step, the medium (20m) and low (60m) resolution bands of Sentinel-2 were sharpened from the nearest high (10m) spatial resolution band.

Class Nomenclature and Reference Data
Concerning the class nomenclature, the aim was to highly analyse land cover of the study area including detailed information for crop types. To this end, we used the nomenclature of CORINE Land Cover (CLC) third level with several modifications for general land cover classes, while geospatial data from Rapid Field Visits (RFVs) of the Greek Paying Agency (OPEKEPE) defined the thematic analysis of arable land classes. An intensive annotation procedure was carried out to produce reference data. Two image interpretation experts manually annotated polygons for 31 different land cover and crop type classes using a variety of datasets including Sentinel-2 images for the year 2016, Google Earth and Bing Satellite very high resolution imaging data, CLC2012 product and crop's geospatial data from Greek Paying Agency. Since representative training samples are one of the most critical components in supervised classification, the experts studied the area thoroughly and noted as many variations of each class as possible (Karakizi et al., 2018b). At the same time, the good practice of keeping sample size per class relative to each class's respective occurrence in the study area was also considered when creating the reference data. However, this practice was impossible to be achieved for arable land classes, since reference data for those -hard to photo interpret-crop classes were limited by the availability of RFV's geospatial data.

Classification Algorithm and Features
A Random Forest (RF) classifier was used for the implementation of the proposed methodology. RF classifier as further analysed in the Introduction section, is considered a robust machine learning classification method that has been successfully applied in recent similar studies (Defourny et al., 2019;Inglada et al., 2017Inglada et al., , 2015Xie et al., 2019). The choice of parameters for the RF classifier has been proven not very sensitive for land cover tasks (Pelletier et al., 2016). For this study we implemented an RF classifier with Python 3.7.0 using scikit-learn and 100 trees. All experiments were executed on a server running Ubuntu 18.04, with an Intel(R) Core(TM) i7-5820K CPU at 3.30GHz and 48 GBs of RAM.
Dataset B was formed by deriving the spectral features from dates/images of Dataset A that fulfilled specific temporal metric conditions. Distinct remote sensing stages, especially for cropareas, have been reported to be defined by key dates linked with spectral features' statistical derivations, like maximum (max) and minimum (min) value of NDVI and Red band (Lambert et al., 2016;Matton et al., 2015). We experimented with multiple combinations of temporal metrics as conditions, like min, max and median of the available spectral bands and the final, most efficient set, comprised of: maxRed, minRed, medianRed, maxNDVI, minNDVI, medianNDVI. Dataset B was formed by extracting the values of the nine spectral features from Dataset A, from the respective date that fulfilled the condition posed by each temporal metric, per pixel. As a result, Dataset B consisted of 54 layers in total.
Dataset C was formed by calculating temporal metrics on the Dataset A time-series and using directly those metrics as input features for the classification. We experimented with multiple combinations of temporal metrics as features and the final, most efficient set, comprised of min, max, median, mean, standard deviation of all nine spectral bands and indices. By this way Dataset C consisted of 45 layers in total.

EXPERIMENTAL RESULTS AND DISCUSSION
In this section a thorough quantitative and qualitative analysis is presented based on the results obtained after applying the same classification framework on Datasets A, B and C.

Accuracy Metrics
All models were trained and validated using exactly the same training and testing areas, in order to allow valid comparison between datasets. The independency of the validation procedure was secured by splitting the reference data at the polygon level, using 65% for training and 35% for testing. An additional validation procedure for Dataset A, was conducted by splitting the reference data randomly at pixel level, holding the same ratio for training and testing, towards quantitatively assessing the impact of the spatial correlation on accuracy metrics' overestimation.
The validation of the land cover classification experiments was quantitatively implemented forming confusion matrices at the pixel level. The standard accuracy metrics of Overall Accuracy (OA), User's and Producer's Accuracy (UA & PA) were calculated. Per class F-measure (F1) scores were also calculated as the combined metric of the harmonic mean between UA and PA. Qualitative evaluation was also performed by throughout intensive observation of the produced land-cover and crop-type maps.

Quantitative Comparative Analysis
In Table 1, the resulting F1 rates are presented per class and dataset, while their average rates and the OA of each experiment are presented in the last two rows.
In particular, lower rates (less than 40%) are marked with double underline and highest rate per class (per row) between the different datasets with bold. Scores on the last two rows of The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2020, 2020 XXIV ISPRS Congress (2020 edition) This contribution has been peer-reviewed. https://doi.org/10.5194/isprs-archives-XLIII-B3-2020-1555-2020 | © Authors 2020. CC BY 4.0 License.
Classification on Dataset A, of 90 spectral features, holds the highest rates (with bold) for average F1 and OA, and for most individual classes' F1 scores, compared to the performances of the other two experiments. Dataset C of 45 temporal features performed, in almost every case, better that Dataset B, of 54 spectrotemporal features. Comparing Dataset A and C, the latter comes with half the dimensionality of the first. However, accuracy rates did not differ much, i.e., less than 3% for average F1 and OA. For the experiments in all three datasets, F1 score for the majority of classes exceeded 70%. For all three datasets, very high F1 rates of over 90% were achieved for classes Greenhouses, Broad-leaved forest, Coniferous forest, Dense sclerophyllous vegetation, Wheat and wetland (MRS) and water classes (WCR, WBD, CWT). Sub-classes of the cereal crop family, name Barley and Oat were the only cases of very low rates of under 40%.
In general, individual classes' performance followed the same pattern, concerning accuracy rates derived from the three different datasets. To better illustrate the latter and at the same time assess the classes' sample size contribution, on the per class accuracy produced by the classification model, Figure 2 is provided and analysed. On the three coloured bars, the PA rates are presented per dataset for every studied class. Producer's Accuracy (PA) is an indicator of the map accuracy from the point of view of the map/model maker (the producer). It expresses the proportion of ground truth pixels that were correctly classified by the classification model. The black line with markers in Figure 2 indicates the proportion of sampling size per class, derived as the ratio between sample size of each class and total sample size for all classes. Sample size is examined since RF, like other supervised machine learning algorithms, tend to promote prediction accuracy of prevailing classes in terms of representativeness which in most cases leads to lower accuracies in less represented classes.
Upon more careful examination of Figure 2 it becomes apparent that the majority of classes, at least 20 out of 31, reached PA scores of over 70% for all 3 datasets. Concerning sample size, the top-five better represented classes, i.e., Coniferous forest, Broad-leaved forest, Cotton, Clover and Maize, achieved in most cases very high PA rates of more than 85%. However PA rates of over 85% were also achieved for very low sample size classes (<1%) like Greenhouses, Marshes and Water courses, suggesting that distinct spectral characteristics of specific classes eliminate the sample size bias. This is indeed highlighted by very high performance (>90%) in both F1 (Table  1) and PA rates for all cover classes related with the presence of water (WCR, WBD, CWT). On the other hand for classes with less distinctive spectral behavior and thus having more competitive classes of similar behavior, like Industrial commercial units and Road and asphalt networks, low sampling size seems to have affected performance against the other manmade classes.
Another important remark from Figure 2 is that the greatest differences in PA scores between datasets A, B and C occurred for crop classes (VNY-FLW). This is also the case for the greatest differences in F1 scores (Table 1). For crop classes other than Vineyards, Tobacco and Other grass fodders, classification on Dataset A ensured highest per class accuracy rates. Sample size seems to have played an important role for those classes and especially the arable ones, that their reference data creation, as explained in section 2.3., was limited by the availability of geospatial data provided by the Greek Paying Agency. In this case, classes Fruit trees, Barley, Oat and Fallow suffered both from low representativeness (<1%) and lower per class accuracy rates. Between those, Barley and Oat presented the lowest per class accuracy rates (in most cases <40%) and this fact can be also attributed to the presence of another crop-class of the same family, namely Wheat that had a significantly bigger sample size than the other two cereal classes. Apart from low sample size, for classes Fruit trees and Fallow lower accuracy rates could be also attributed to the fact that they include different covers or crop sub-types within the same class/category.

Omission Errors per Class
In order to have a view of classification mixings between classes, we hereby present and discuss Table 2. Table 2 was calculated by dividing the produced confusion matrix for the most successful experiment, that of Dataset A, with the sum of reference testing pixels (ground truth) for each row. As a result, the table should be read per row (reference labels), so as the diagonal of the matrix expresses the PA percentage of each row/class, and all other cells the omission percentage for that class from classification errors towards all other column-classes (produced labels). Black line borders have been designed in order to visually group classes of the same general land cover family.
Artificial classes that consist of man-made surfaces achieved, in overall, medium to very high PA rates (49-98%) and in most cases they were confused between each other. Photovoltaic units and Greenhouses, characterised by materials of distinct spectral behaviour, held high rates for PA (>80%) and very low omission errors accordingly. The two classes of urban fabric (DUF, SUF) presented mixings of about 15% between each other, as expected. Omission errors towards those two classes, adding up to approx. 30%, presented the class Road /Asphalt networks with a PA of 49%. The artificial terrain class Mineral extraction sites, which is the only non-sealed man-made cover category in our nomenclature, achieved a PA rate of 58%, while presenting a high omission error rate (34%) to class Beaches, dunes, sands, of similar cover materials. Industrial commercial units class had omission errors towards both artificial (DUF, SUF) and bareland (BDS, RCK) classes of high brightness, presenting a PA rate of 60%.
For the 12 classes of the agriculture family PA rates presented a large range (31%-97%). Permanent crops (VNY, OLG, FRT) were classified with PA rates of about 70%-75% presenting omission errors mainly to other crop classes but also natural vegetation classes and especially towards Sparse sclerophyllous vegetation. Arable land classes Wheat, Maize, Cotton, Tobacco and Clover all held very high PA rates (>89%) and presented insignificant omission errors mainly to other crop classes.
Barley and Oat classes of the cereal family, presented lower PA rates, since 36% and 25% respectively, of their reference data pixels used for testing, were classified as Wheat. As also discussed in the previous section, cereal class Wheat had a significantly bigger sample size than the other two cereals. Additionally Other grass fodders class also lost 49% of its ground truth pixels to Wheat. Fallow class, not corresponding to a single crop category, presented lower PA rates and omission errors mainly to Wheat (14%), Cotton (6%), Natural grasslands (9%) and Sparse urban fabric (11%). The latter misclassification case can be attributed to similar spectral behaviour at this scale (10m) since Fallow and Sparse urban fabric classes include both vegetated and non-vegetated covers in an irregular way.
Regarding bareland classes (SVA, RCK, BDS) very few omission cases were recorded and so PA rates exceeded 85%. Unsurprisingly Beaches, dunes, sands class presented omission errors of 4% to Water courses class, since its reference data also included sand along river banks. Bare rocks class had a 7% of its reference testing pixels classified as Mineral extraction sites class, strengthening the spectral relationship, mentioned in previous paragraphs, between those different-family but similarmaterial classes.
Wetland (MRS) and water classes' (WCR, WBD, CWT) reference testing pixels were classified with very high accuracy rates of over 91%, highlighting the spectral distinctness of water-related covers especially when exploiting the infrared spectrum.

Overestimating Accuracy with Spatial Dependency on Training and Testing Sets
Machine learning and remote sensing best practices for validating classification outputs suggest splitting testing and training datasets in a manner that ensures the spatial independency of the two sets. This is also the case for land cover and crop type mapping applications (Defourny et al., 2019;Inglada et al., 2017;Stoian et al., 2019). This splitting strategy was followed also in this paper for the core experiments presented and analysed in the previous sub-sections.
However, in order to quantitatively assess the over-estimation effect on accuracy metrics, when training and testing data are not spatially independent, an additional experiment was carried out on Dataset A by splitting the reference data at pixel level instead of polygon level.  Table 3 indicates that all metrics exceeded the remarkably high rate of 95% when splitting at pixel level presenting significant raises compared to the main polygon-split experiment. In particular, OA for the complementary experiment exceeded 98% presenting an increase of more than 5%, while average per class metrics of PA, UA and F1 presented rises of over 10%. Rises were quite expected since splitting at pixel level involves considering pixels of the same polygon for both training and testing, which would probably be very similar spectrally.
After all, producing land cover maps covering large extents, such as those of one satellite scene, with not-automatically produced reference data from other LC products, is usually characterized by a small proportion of reference data compared to the whole area to be classified. In such a way the impact of spatial correlation when validating, is not to be ignored in order to provide the users an accuracy estimation that is closer to reality for the whole map.

Map Validation and Discussing Challenging-to-map Areas
The final step of the proposed methodology is map production. Three maps were produced, one for each dataset, applying the same RF framework. For predicting the map products the whole reference dataset was utilized when training the RF classifier, aiming at maximum exploitation of all available information. The general visual inspection on the produced maps did not reveal coarse differences between the different products. The most successful experiment i.e., the one based on Dataset A, is presented in Figure 3 (left hand side). Additionally, zoomed-in areas from all three maps are presented for visual comparison in Figure 3 (right hand side). This allows the direct comparison of the performance especially in challenging-to-map areas.
On the eastern part of the map (Figure 3, left) one can observe that the urban and suburban areas, i.e., cities of Trikala and Karditsa as well as many villages scattered across the Thessalian plain have been correctly detected and mapped (indicative classes: DUF, SUF, ICU) by the proposed classification framework. Forest classes (BLF, CNF) were correctly reported on upland and mountainous areas, while Mediterranean-type evergreen sclerophyllous vegetation classes (DSV, SSV) were mapped on lower terrains.
Regarding the qualitative evaluation over agricultural areas, the Thessalian plain was expectedly classified as various crop classes. A vast dominance of classes Cotton, Wheat and Maize cultivations is evident on the map which agrees with the agricultural information provided in subsection 2.1. Additionally, large inland water areas such as the man-made lakes of the study area, the salty-water lagoon of Messolonghi and its surrounding wetlands have been successfully detected as the respective classes (WBD, CWT, MRS).
Higher altitude regions of mountain massifs, mainly represented by the spur of Pindus, seem in most cases correctly illustrated with the presence of Natural grasslands and bareland classes (SVA, RCK). However a closer look on bareland areas (whitish colors) even at this scale, reveals certain omission errors towards man-made classes (redish/pink colors) that are further illustrated on the zoomed-in parts in the right side of Figure 3.
In particular the example of mount Lakmos (alt. 2100m) of the Pindus range is presented here. This area consists of bareland covers and usually snow cover across the year. Therefore it is depicted in the imagery data with high brightness values across the spectrum, which is also the case for many other alpine regions.
By examining the zoomed-in parts from the maps of Dataset A, B and C, classification errors are apparent for all. The zoomedin map of Dataset A present the fewest commission errors toother than bareland classes-among the three maps. However, there are still considerable regions covered with urban and suburban classes (DUF, SUF) but also class Cotton of the crop family. The artificial classes' mixings with bareland can be attributed to high brightness values for both types. As far as class Cotton is concerned, its cultivation practices require high humidity levels on the soil, so it is possible that alpine vegetated regions covered by moss and lichen were confused with this class. More classification errors can be observed for Dataset B and C. Larger bareland areas are erroneously classified to even more in number artificial classes, namely Dense urban fabric, Sparse urban fabric, Mineral extraction sites and Greenhouses. These remarks are in accordance with per class accuracies, since classification in Dataset A resulted in Sparsely vegetated areas PA rates of 86%, 85%, 80% and that of Bare rocks of 90%, 88%, 87% for datasets A, B and C respectively.
In overall, visual inspection on the produced three maps agrees with the high classification OA rates (>90%) achieved by all datasets. Additionally, analysis on the challenging-to-map example, indicated also qualitatively Dataset A as the most successfully classified dataset.

Classification map on Dataset A
Zoom-in areas Figure 3. Classification map of the study area from experiment on Dataset A (left) and zoom-in region of mount Lakmos from the maps of all three datasets along with Bing Satellite imagery of the same region (right)

CONCLUSIONS
In this paper we assessed the contribution of spectral and temporal features for detailed land cover and crop type mapping using annual Sentinel-2 data for the year 2016. Three datasets were created and benchmarked using the same RF classification framework. The first baseline dataset (Dataset A) consisted of 90 spectral features, while the other two were formed with spectrotemporal (Dataset B) and temporal features (Dataset C), reducing the initial datacube dimensionality by 40% and 50% respectively. Quantitative evaluation demonstrated the efficiency of the proposed methodology by achieving results of high overall accuracy of over 90% in all experiments. Comparative quantitative and qualitative validation highlighted that in most cases Dataset A yielded modestly better results than the other two datasets that employed temporal metrics.
The choice of using temporal metrics for this work has been encouraged by recent similar studies and especially crop mapping ones, with the aim of allowing independency from the crop calendar and dealing with the cropland diversity across landscapes (Defourny et al., 2019;Egorov et al., 2018;Waldner et al.,2017). In this work, the study area consisted of one Sentinel-2 tile and thus cropland diversity and agro-climatic gradient did not ranged significantly. To this extent, the classification experiment based only on temporal metrics (Dataset C) presented lower rates by a 2-2,5% for OA and average F1, compared to the initial experiment on full-set spectral features of Dataset A. Nevertheless relative low differences in accuracy for experiments performed, encourage the use of temporal metrics for similar mapping tasks and especially for large geographical areas, contributing to reduced data dimensionality and thus less demands on hardware resources.
Further analysis on per class classification errors indicated that sample size along with the distinct spectral characteristics of each class had a significant impact in accuracy rates. Especially for arable crop classes, that present high spectral similarities among each other, reference data availability emerges as a crucial prerequisite for accurate classification. Towards this end, the free/open access to geospatial reference datasets like the annual declarations on the Land Parcel Identification System (LPIS, for EU) should significantly contribute to accurate and detailed land cover mapping at an annual basis.
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLIII-B3-2020, 2020