ARCHAEOLOGICAL PROSPECTION WITH CORONA AND WV-3 SATELLITE IMAGERY OF THE ARCHAEOLOGICAL SITE OF ZAR TEPE (UZBEKISTAN)

: The use of different data from satellite platforms for archaeological prospecting and remote sensing has been applied since the end of the 20th century. Although the current use of drones with different visible and multispectral sensors for small areas has partially replaced in some cases the use of this type of satellite information due to its higher spatial resolution. The historical importance of satellite imagery is essential to find out about and compare the transformations of the archaeological landscape of the last 60 years since the CORONA satellite program started in 1960. In this paper we propose the evaluation of a proposal for the automation of processes of two photographic reconnaissance correlative satellite programs CORONA (1960-1972), HEXAGON (1971-1986) declassified since 1995 and 2011 respectively, and the commercial satellite WorldView-3 (WV3) (2014) for use in the detection of buried archaeological structures at the archaeological site of Zar Tepe in the southeast of Uzbekistan. This is a site located in the Surkhan Darya region very little known between the first century BC and the fourth century AD. This methodology is part of the IPAEB project (International Pluridisciplinary Archaeological Expedition to Bactria) led by the University of Barcelona, the University of Salamanca and recently in 2019 the University of Zaragoza. IPAEB was started in 2006 in the South of Uzbekistan and is currently trying to explore the urban planning of the Zar Tepe archaeological site and the elements that make up its natural physical environment: evidence of communication routes, smallholdings, irrigation channels, fences and sources of raw materials.


INTRODUCTION
Recent improvements in the fields of remote sensing and computational processing have made possible the development of new methodologies applied to the use of declassified program data such as CORONA (1960)(1961)(1962)(1963)(1964)(1965)(1966)(1967)(1968)(1969)(1970)(1971)(1972) and HEXAGON (1971HEXAGON ( -1986 photographic reconnaissance satellites to document archaeological sites (Hammer et al., 2019(Hammer et al., , 2022Casana, 2020). This type of information has been easily accessible since its declassification in 1995 (CORONA) and 2011 (HEXAGON) and made available to the public through the US Geological Survey. These satellite reconnaissance programs complement each other in their global coverage and resolution and are even able to extract DEMs thanks, for example, to the stereoscopic pairs of the CORONA KH-4. These types of features are very useful in archaeological landscapes that have changed in the last 50 years due to anthropogenic or natural factors (Angas, 2019(Angas, , 2021. Likewise, satellites such as WorldView-3 (WV3) (2014) and satellite constellation such as Pléiades Neo carry very highresolution optical instruments (30 cm) which could be useful for identifying subsurface archaeological structures (El-Behaedi, 2021). The evaluation and comparison of both data sources to correct and process remote sensing data programmatically allow the creation of reproducible workflows to obtain and analyse vegetation cropmarks of archaeological areas. Research in this area of Uzbekistan began in 2006 with the International Pluridisciplinary Archaeological Expedition to Bactria (IPAEB) group from the University of Salamanca, the University of Barcelona and recently in 2019 in the University of The objective of this paper is to present the evaluation of the methodology applied to the automation of processes and detection of buried urban structures in the archaeological site of Zar Tepe, in the southeast of Uzbekistan to provide more information about this site given the scarce archaeological documentation compared to other settlements in the Surkhan Darya Valley.
The resulting cartography will serve as support in the field prospecting campaigns of 2023 and 2024.

STUDY AREA
The documentation focuses on the archaeological settlement of Zar Tepe (Golden Hill), in southeastern Uzbekistan (Figure 1). This site is located in the Surkhan Darya region, 26 kilometres north of Termez. Its foundation is dated with little precision between the first century BC and was occupied until the fourth century actively, and until the sixth century AD sporadically (Zavyalov, 2008). It is located on the east bank of the Amu Darya River, close to the Surkhan Darya mountains. It is part of the Afghan-Tajik plain, composed of metamorphic rocks from the Precambrian-Paleozoic era and sedimentary rocks from the Triassic. The nearest city is Termez, one of the most important in northern Bactria/Tokharistan, a historical region that includes present-day Tajikistan, southern Uzbekistan, and northern Afghanistan (Sánchez del Corral et al. 2002). The city was built at the end of an irrigation canal that stems from the Surkhan River, in a strategic position on the route that connected Bactria and Sogdiana (Termez with Samarkand). Zar Tepe illustrates the different stages of evolution of a Bactrian city in the Kushan period and has been the subject of archaeological excavations since the mid-20th century (Zavyalov, 2008) but most of this archaeological site has not been investigated; the higher levels have been prioritised but there are important questions related to the founding date. Its evolution and historical meaning are still uncertain and need to be studied in order to gain a broader and more accurate understanding of Zar Tepe and the Bactrian cities in the Kushan and post-Kushan times.
The city has a square plan 400m on each side, surrounded by a wall provided with an interior corridor, semicircular towers 6m in diameter, spaced every 34m and a moat. At its N-NE angle, a square fortification 120m on each side was built on an artificial platform and isolated by a citadel moat. The most important constructions were concentrated around the main road that crosses the city from East to West, and other adjacent roads; the SE sector housed the residential neighborhood, with houses of different sizes (between two and eight rooms) and an isolated enclosure interpreted as a military fort; other civil buildings and craft workshops were located to the northwest. Although pottery kilns are not known at Zar Tepe, the documented ceramic contexts are notable and respond to a varied repertoire of forms in all functional categories (kitchen, storage, tableware, common ceramics, oil lamps), which shows certain analogies with the ceramics of Termez. Zar Tepe's function as an important economic, religious and political center in the Surkhan Darya valley between the first and third centuries is evidenced by the residential palace from the post-Kushan period, comprised of two columned rooms (the Governor's audience rooms) provided with a perimeter bench and another in a central position corresponding to that of the throne. The Buddhist cult is represented by a sanctuary located in the north-central sector, which provided fragments of a golden Buddha statue and a bronze vessel for the sacred fire.
The objectives proposed in this research in the South of Uzbekistan and especially in the Surkhan Darya valley are: 1) Investigate the urbanism of Zar Tepe and the elements that make up its physical-natural environment (evidence of communication routes, smallholdings, irrigation channels, fences, sources of raw materials). In this sense, this investigation has been approached through a macro-spatial study using noninvasive remote sensing techniques and verification by superficial prospecting in the following archaeological campaigns.
2) Reassess the occupational sequence of Zar Tepe to confirm or refute the chronologies and historical evolution attributed up to now, verify whether its foundation was promoted by the Yuezhi chiefs after their settlement north of the Oxus and at what time, or by the Kushan kings in the 1st century AD. It will be approached through the prospecting and excavation of a sector of the citadel, scarcely explored by earlier research, and which could house the oldest remains. The extent of the excavation in 2023 and 2024 will depend on the results of the remote sensing survey.
3) Explore the existence of possible potteries, verify whether the anomalies shown by the satellite images in the extreme Northwest (artisanal area) correspond to ceramic kilns. This will be addressed through remote sensing and superficial verification prospecting. If ovens are located, at least one will be excavated and documented, to which the photogrammetric record and Carbon 14 analysis will be applied.
4) Document, contextualise and formally and stylistically categorise the ceramic products found in the excavations; compare this data with the absolute dates and the chronotypological series created by the IPAEB for the contemporary ceramics from Termez to check if the results can be extrapolated to other Bactrian settlements.

Satellite images
Satellite imagery such as CORONA (Angás 2019(Angás , 2021, Hexagon and WorldView 3 have been used for this research. CORONA images are derived from KH (Key Hole) 1 -4 panoramic cameras. These were developed within US reconnaissance satellite programmes (1960 -1972) to collect lower resolution photographic coverage at a broad scale (Hammer, FitzPatrick and Ur, 2022). The best spatial resolutions for each satellite ranged from 12 to 1.8 metres. Lenses from these cameras varied from 0.5 to 0.7 µm, detecting the visible fraction of the electromagnetic spectrum, with 0.62 and 0.65 µm as optimum (Dashora, Lohani and Malik, 2007).
Hexagon images were captured with KH-9 cameras, taking pictures between 1971 and 1984. The camera specs were the same as CORONA KH-4A and KH-4B cameras, the last of the CORONA series. The best spatial resolution of KH-9 cameras ranged from 0.6 to 1.2 metres.
Initially, the CORONA and Hexagon imagery was used, between 1969 and 1978 (table 1); these images have a single band and resolutions below 2m. They correspond to a period with less anthropisation, where the appearance of soil marks depends on atmospheric and light conditions. These images were acquired using the USGS EarthExplorer browser (National Geospatial-Intelligence Agency -NGA) divided by frames which are sorted by letters. Due to errors in the georeferencing of the images in Declass 1 (1996), they were filtered with an area of interest (AOI) of 100,000km 2 , thus ensuring the selection of all available images. Subsequently, the WorldView-3 imagery was used, which was acquired with eight bands between the blue and near-infrared spectral regions. These have a spatial resolution of 0.31cm in the panchromatic band in NADIR and 1.24 m in the multispectral bands. Three images have been acquired that correspond to three outstanding moments of the phenological cycle of cereal crops, i.e., autumn/birth, spring/greening and summer/senescence.

Image processing
The second phase of the methodology consists of the treatment of the images. On one hand, CORONA scenes must be individually georeferenced because they do not have coordinate reference systems. Additionally, they contain large spatial distortions which are more important towards the perimeter. QGIS (QGIS Development Team, 2022) has been used together with the Georeferencer plugin to georeference each image.
Initially, the image is clipped to the AOI of the archaeological site using the GDAL warp tool (GDAL, 2023) to reduce the area where ground control points are located.
Using expert criteria, stable control points are placed between the CORONA image and the panchromatic band of the WorldView-3 image (e.g. crossroads, irrigation canals, anthropic constructions, etc.). Newly referenced points are then obtained in the EPSG:32642 WGS 84 / UTM zone coordinate system 42N.

Figure 2. Workflow scheme.
On the other hand, the WorldView-3 images are fixed programmatically using a Python script which performs the following process for each of the images using its metadata ( Figure 2): (1) transformation of the TIL image to GeoTIFF format; (2) application of the conversion from relative radiometrically corrected image pixels (or digital numbers, DN) to radiance at Top Of Atmosphere (TOA) using the Absolute Radiometric Correction (ARC) transformation; (3) correction of the image atmospherically with the Second Simulation of the Satellite Signal in the Solar Spectrum (6S) (Vermote et al., 1997). The output for each image is a corrected panchromatic and multispectral WV3 image.
The first step is performed to compress the image, and the remaining steps for applying the atmospheric correction to obtain reflectance at the Bottom of Atmosphere (BOA). The ARC transformation formula (Equation 1) is derived from the WV3 technical sheet (add reference), where X is the value from the original pixel; GAIN and OFFSET are the variability in electronic gain and offset parameters derived by band to fix non-uniformity inside the image; abscalfactor and effectiveBandwidth are flightdependent parameters included in the IMD metadata file. The script iterates through all metadata files to retrieve the required data in the process by using the WV3 naming convention.
To complete the atmospheric correction, the TOA radiance values are transformed to BOA reflectance by subtracting the atmospheric effect. The 6S model is applied with both the Py6S python module (Wilson, 2013) and the Google Earth Engine Python API to retrieve the required atmospheric parameters in the radiative transfer equation (2): where L is TOA radiance; Lp is the path radiance or the amount of energy captured by the sensor from the atmosphere interaction; τ is transmissivity (absorption transmissivity divided by the scattering transmissivity), Edir is the direct solar irradiance and The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-M-2-2023 29th CIPA Symposium "Documenting, Understanding, Preserving Cultural Heritage: Humanities and Digital Technologies for Shaping the Future", 25-30 June 2023, Florence, Italy Edif is the diffuse solar irradiance. All the prior terms except L value are solved by the 6S model. Once the images are corrected, the pansharpening process is performed between the WV3 corrected panchromatic band and the multispectral corrected image using the weighted Brovey technique, which is applied with the GDAL pansharpen tool (GDAL/OGR, 2022) using the bands and percentages defined by Belfiore, 2016. The HPF Resolution Merge (Gangkofner, Pradhan and Holcomb, 2008) pansharpen method was also applied ( figure 3). However, the results do not improve the first technique in terms of in terms of colour and sharpness (Zhang and Mishra, 2012). Image analysis was performed by professional visual quality evaluation, due to the lack of quantitative methods which are able to provide an accurate assessment (Zhang, 2008).

Create web map tool
To ensure all the professionals involved in the project could analyse all the data, a web visor has been created from scratch with html and JavaScript code. GIS software programs are the best option to represent and combine spatial data and georeferenced images, but these are difficult to manage. The proposed web visor has a simple interface with a browse panel where the user can select/unselect the images to view on the web. The spatial data is represented in a web map programming with JavaScript Leaflet library, and there is an option to draw the cropmarks on the screen and download it in geoJSON format.

Create neochannels
Once processed, both the World-View-3 and CORONA images have been subjected to a series of band/multiband operations using a second Python script. When the operation requires the WV3 second infrared band (B8), the original image will be used. This band could not be included in the pansharpen operation because it is outside the spectral range of the panchromatic band. The first multiband operation is the calculation of twenty-two spectral indices (Table 3), combining the bands that allow the highlighting of the differences in moisture, soil components and plant vigour.
The second operation consists of calculating the principal component analysis, as many times as band combinations can be obtained without repetition (a minimum of three and a maximum equal to the number of bands in the image). From each combination, the three components with the highest eigenvalues are obtained, i.e., with the largest variations around their axes.
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-M-2-2023 29th CIPA Symposium "Documenting, Understanding, Preserving Cultural Heritage: Humanities and Digital Technologies for Shaping the Future", 25-30 June 2023, Florence, Italy

Detect crop marks
The above products (figure 4) are analysed together with the web mapping tool created from scratch to detect cropmarks on the field. The results obtained will allow a verification by means of a topographic stakeout in the field in the determined areas. In conclusion, CORONA imagery is usually seen clearly when the sun is low in the sky and shadows accentuate small changes in the terrain . Regarding the WorldView-3 imagery, the selection of images in the three phenological cycles has made it possible to determine that this type of crop is the one that best reflects the underground structures due to its root depth, with greater exposure to water stress produced by the presence of elements in the subsoil.

CONCLUSIONS
After this preliminary research study, we cannot affirm that any post-processing works better than another for the visualization the reconstruction of buried archaeological structures, but that all The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVIII-M-2-2023 29th CIPA Symposium "Documenting, Understanding, Preserving Cultural Heritage: Humanities and Digital Technologies for Shaping the Future", 25-30 June 2023, Florence, Italy of them complement each other, allowing different soil marks to be observed in each one of them. It is clear that, in this specific case, the best results have been obtained with the WV3 satellite and its multispectral images. Its spatial resolution, multitemporality and its spectral resolution allow us to obtain different products that improve the visualization and interpretation of the remains. Thanks to the vegetation indices we can better interpret the wall and the wet areas or areas with stagnant water. According to the PCAs, the city's road network becomes clearer, a fact that does not happen with the indices. If in the case of Termez (Angas et al. 2019(Angas et al. , 2021 we obtained excellent results with the CORONA program nevertheless the same does not happen in the case of Zar Tepe. The poor image sharpness comes as a result of the lack of spatial resolution of the images from the satellite missions. Consequently, only the 1973 image offers a higher spatial resolution (0.6-1.2 m) than the rest (1.8 m). This fact allows us to see more clearly the inner wall of the citadel and its configuration. Although the post-processing of the images is automated, the work of reviewing the data with more than 200 PCAs requires significant verification time. Subsequently, in the future we will try to apply artificial intelligence at least for the automatic recognition of certain anthropic elements. Regarding the workflow presented, it is a suite of programs to process and create a number of derived products which are able to detect cropmarks on the field. In addition, the custom-created web mapping is a useful tool to allow professionals with no GIS knowledge analysis to analyse all derived products and search crop marks remotely. A future development could be added a new product derived with new neuronal network clustering techniques (Conner, Carter and Blackadar, 2021).