TOWARD AUTOMATIC GEOREFERENCING OF ARCHIVAL AERIAL PHOTOGRAMMETRIC SURVEYS

Images from archival aerial photogrammetric surveys are a unique and relatively unexplored means to chronicle 3D land-cover changes over the past 100 years. They provide a relatively dense temporal sampling of the territories with very high spatial resolution. Such time series image analysis is a mandatory baseline for a large variety of long-term environmental monitoring studies. The current bottleneck for accurate comparison between epochs is their fine georeferencing step. No fully automatic method has been proposed yet and existing studies are rather limited in terms of area and number of dates. State-of-the art shows that the major challenge is the identification of ground references: cartographic coordinates and their position in the archival images. This task is manually performed, and extremely time-consuming. This paper proposes to use a photogrammetric approach, and states that the 3D information that can be computed is the key to full automation. Its original idea lies in a 2-step approach: (i) the computation of a coarse absolute image orientation; (ii) the use of the coarse Digital Surface Model (DSM) information for automatic absolute image orientation. It only relies on a recent orthoimage+DSM, used as master reference for all epochs. The coarse orthoimage, compared with such a reference, allows the identification of dense ground references and the coarse DSM provides their position in the archival images. Results on two areas and 5 dates show that this method is compatible with long and dense archival aerial image series. Satisfactory planimetric and altimetric accuracies are reported, with variations depending on the ground sampling distance of the images and the location of the Ground Control Points.


INTRODUCTION
Archival aerial photogrammetric surveys were initially acquired by mapping, cadastral or military agencies for topographic map generation.Stereoscopic configurations were adopted so as to offer 3D plotting capacities.These surveys have been a common practice in many countries over the last century.Recently, several countries have digitised their film-based photos (e.g., >3 millions images in France), and facilitated their access through spatial data infrastructures and web services with basic metadata and visualisation capacities.These images are a unique yet unexplored means for long-term environmental monitoring and change analysis, while they chronicle Earth surface evolution in a comprehensive way.Their use in remote sensing workflows requires a fine georeferencing step between dates so as to generate orthoimages that can be spatially compared.In addition, when a photogrammetric method is considered (dense matching following interior and absolute raw image orientation), Digital Surfaces Models (DSMs) can even be obtained.
Until now, satellite images and archival topographic maps have been used to characterise land-cover dynamics.On the one hand, satellite images, such as Landsat and SPOT sensor families, are compatible with automatic analysis requirements: their radiometry and spatial resolutions are stable over time.For instance, Landsat images have been widely used in order to monitor urban land-cover change (Song et al., 2016).However, the types of changes that can be described are very limited by their low spatial resolution (30-60 m), the length of the time series (30-40 years) and to 2D (Hermosilla et al., 2015).On the other hand, topographic maps are available in many countries since the 18 th century, but they are difficult to analyse automatically.Thus, approaches that characterise land-cover evolution are heavily based on Volunteered Geographic Information (Perret et al., 2015) and their temporal resolution is very coarse: the objects that can be extracted (buildings, roads) are difficult to date precisely (Leyk et al., 2006, Herrault et al., 2013).
In contrast, archival aerial images are perfectly tailored to characterise long-term land-cover changes (Giordano et al., 2017).The time series are (i) long (first images: 1920s), (ii) relatively dense temporally (every 2-5 years), (iii) acquired at very high spatial resolution (<1-2 m) and, (iv) more remarkably, the survey conditions enable to generate 3D information.However, such advantages are accompanied with very particular characteristics.They are usually composed of an unique panchromatic channel, and their radiometry can be very variable.The digital images are generated from old analog photography, acquired with different cameras for which only coarse image localisation was stored.In addition, artefacts and noise are often observed: radiometric issues heavily impact the geometric accuracy of the orthoimages (Micheletti et al., 2015), and subsequently the computation of the DSMs.Moreover, analysing these time series is rather challenging since images exhibit highly heterogeneous spatial resolutions, with very different acquisition conditions.Useful metadata is generally distributed, such as length and width of the analog photos, focal length of the camera and more rarely camera calibration reports.Flight plans of the aerial surveys give very coarse image localisation (∼100 m) and are often distributed with the images.
Very few studies have proposed dense long-term environmental monitoring and change analysis based on archival aerial images.The challenges can either come from the lack of fine georeferencing methods compatible with the large amount of data or the difficulty to analyse such heterogeneous time series.This paper focuses on the georeferencing stage and proposes a coarse-to-fine solution for that purpose, enabling both to provide comparable orthoimages and DSMs.State-of-the-art is provided in Section 2, while Section 3 describes the proposed approach.Results are presented in Section 4, and conclusions are drawn in Section 5.

FINE GEOREFERENCING OF ARCHIVAL IMAGES
In this paper, the hypothesis is made that a very coarse image location is supplied for each archival aerial image (∼100 m).It therefore focuses on fine georeferencing of such data into a cartographic system.A precise accuracy is expected so as to be able to perform comparison between the products that can be derived (radiometry, height, labels).This paper proposes to categorise the methods in two approaches: non-photogrammetric and photogrammetric approaches.Non-photogrammetric approaches find an empirical transform that gives the best match between images: most of the standard registration methods fall into this category (Zambanini and Sablatnig, 2017).Conversely, photogrammetric approaches aim at modelling the geometry of all images.A DSM can thus be computed.
An overview of existing approaches according to their spatial and temporal characteristics is provided in Figure 1.First, one can notice that most of the studies are based on photogrammetric approaches (6 non-photogrammetric vs. 16 photogrammetric).None of them proposes an automatic solution.For photogrammetric approaches, a majority of studies concern less than 5 dates and 100 km 2 .By increasing size, these papers are: (Meyer et al., 1996, Nebiker et al., 2014, Walstra et al., 2004, Mondino and Chiabrandoa, 2008, Véga and St-Onge, 2008, Kadmon and Harari-Kremer, 1999, Fox and Cziferszky, 2008, Cardenal et al., 2006, Necsoiu et al., 2013, Aguilar et al., 2013, Ayoub et al., 2009).(Micheletti et al., 2015) (9 dates), (Nurminen et al., 2015) (10) and (Korpela, 2006) (28) are the only authors that have processed more than 5 epochs.(Nurminen et al., 2015) consider relatively larger area (∼ 600 km 2 ) for DSM generation in forested areas.Same conclusions (restricted epochs and area size) can be made with non-photogrammetric approaches (Ellis et al., 2006, Jao et al., 2014, Chen et al., 2016), although some of them try to process larger areas (Asner et al., 2003, Galster et al., 2008, Ford, 2013).The time interval and image age are not limiting factor: the average time interval is 30 years and the starting year around 1945-1955. Spatial (study area) . Spatial (study area) and temporal (number of epochs) lim-itations can be due to many factors and probably mainly to the significant manual implication in the fine georeferencing process.As this paper focuses on fine georeferencing, one can try assessing the amount of manual labour.The non-photogrammetric and photogrammetric approaches are evaluated separately (Sections 2.1 and 2.2, respectively).

Non-photogrammetric approaches
Fine non-photogrammetric georeferencing methods are usually based on the application of a polynomial transform to each archival image of the aerial survey.These transforms are estimated using necessary ground references that can be obtained with ancillary georeferenced orthoimages or GPS field surveys.This ground reference is usually a point, namely Ground Control Points (GCPs).These approaches imply the independent identification of dense ground reference for each archival image of the survey: 20 GCPs per archival aerial image for (Ellis et al., 2006), between 10 and 47 GCPs for each image for (Ford, 2013), 56 GCPs for (Asner et al., 2003), and, at least 6 GCPs per image for (Galster et al., 2008).The manual acquisition of the ground coordinates of the GCPs and the estimation of the transform are done using GIS softwares.The ground reference is acquired using recent satellite imagery (Ellis et al., 2006, Ford, 2013), recent national aerial orthoimages (Ellis et al., 2006, Galster et al., 2008) or GPS survey on road intersections (Asner et al., 2003).The registration model can be for example a first order affine transform (Galster et al., 2008) or a second order polynomial transform (Ford, 2013).The planimetric accuracy obtained is usually found compatible with multi-temporal analysis (<1 m, (Asner et al., 2003, Galster et al., 2008) -<2 m (Ford, 2013)).
Finding corresponding GCPs between recent and historical datasets is a complex task, in particular in areas facing significant changes.(Chen et al., 2016) argue this issue makes the process of all the available archival aerial surveys impossible.As a consequence some authors have proposed to use standard automatic tiepoint extraction and matching techniques (Lowe, 1999).However, the tie-points are only extracted between the archival images themselves to reduce the number of manual GCPs needed to estimate the polynomial transform.No tie-points have been extracted between the archival images and a ground reference (satellite or aerial orthoimage).(Jao et al., 2014) propose to estimate each polynomial transform (6 affine transform parameters) using tiepoints between archival images and a more restricted numbers of manual GCPs (6 GCPs for 6 archival images).However, the planimetric accuracy obtained is not suitable (17.6 m).Using the same approach that combines automatic archival tie-points and manual GCPs (Chen et al., 2016) obtain a 1 m planimetric accuracy using a 2D transform (20 manual GCPs for 12 images).
Non-photogrammetric methods are suitable to obtain georeferenced archival orthorectified images with a suitable planimetric accuracy.However, even if some improvement have been proposed, compulsory manual interaction remains a major limitation to move to larger spatial or temporal scales.The extraction of tie-points between archival images and reference orthoimages should be investigated.The lack of initial coarse georeferencing for the archival aerial images and the possible important landcover changes are the main challenges.First attempts have only been tackled for approximate georeferencing and under limited land-cover changes (Kim and Pollefeys, 2008).

Photogrammetric methods
Several photogrammetric softwares are adapted to process archival aerial image surveys (Pierrot Deseilligny andCléry, 2011, Nebiker et al., 2014).A typical workflow that generates true orthoimages and DSMs is given in Figure 2. In the following, the current automation level is assessed for each step.Interior orientation The analog photographs were deformed and thereafter scanned.The interior orientation consists in establishing a mathematical relationship for each image between the camera frame coordinate and the image coordinate system (analog photo).This is done using the fiducial marks of each image.Almost all approaches identified in the literature are applied to aerial photographs with existing fiducial marks.Calibration certificates can be used to retrieve the fiducial mark physical coordinates (Nurminen et al., 2015, Korpela, 2006, Micheletti et al., 2015, Véga and St-Onge, 2008), but they are barely available today.This is a first obstacle to move to larger temporal intervals.As a result, some authors propose to use coarse calibration information (calibration models).(Aguilar et al., 2013) only use the focal length information usually saved or mentioned on the photographs, whereas (Nocerino et al., 2012) absolve themselves from calibration certificates by recovering an approximate information from the images.The authors do not mention how to retrieve image coordinates (pixels) of fiducial marks.This could be done manually.However, some photogrammetric softwares (e.g., MicMac (Pierrot Deseilligny and Cléry, 2011)) propose semi-automatic fiducial mark detection methods based on pattern matching.The manual digitization of the fiducial marks is only necessary on a single archival image of the aerial survey.

Relative image orientation
The relative image orientation is carried out using automatically detected tie-points and the interior orientation.Manual intervention is not necessary at this stage, tie-point extraction and matching being a solved problem in photogrammetric computer vision.All the authors use automatic tiepoints extraction and matching.During this stage, the calibration of the camera of an aerial survey is usually questioned and reestimated.This self-calibration is done automatically using the tie-points (Nurminen et al., 2015, Cardenal et al., 2006, Necsoiu et al., 2013, Fox and Cziferszky, 2008, Walstra et al., 2004).Nevertheless, it must be kept in mind that results of self-calibration often remains an approximation because the poorly geometrically constrained configuration of the image sets.The relative image orientations are estimated using the tie-points and the camera calibration.The relative image orientation has a great impact on the planimetric and altimetric accuracies of the orthoimages and the subsequent DSM.Therefore, some authors chose to add manual intervention during the relative image orientation stage: tie-points are for instance manually extracted in (Micheletti et al., 2015).Multi-temporal tie-points can even be inserted (Cardenal et al., 2006).

Absolute image orientation
Producing an absolute image orientation in a cartographic system requires coordinates of ground references.The most common method is to measure Ground Control Points (GCPs) and to determine their image coordinates in the relevant aerial photographs.In the literature, GCPs are measured manually using ancillary orthoimage and a Digital Terrain Model (Nurminen et al., 2015, Nocerino et al., 2012, Necsoiu et al., 2013, Véga and St-Onge, 2008, Ellis et al., 2006, Fox and Cziferszky, 2008), field GPS (Micheletti et al., 2015, Cardenal et al., 2006, Walstra et al., 2004), or with both techniques (Aguilar et al., 2013).The number of GCPs needed is important compared to the size of the study sites: around 45 GCPs for less than 10 archival aerial images in (Nurminen et al., 2015) or (Aguilar et al., 2013), or 169 GCPs for a 20 km 2 area in (Micheletti et al., 2015), 10 to 47 GCPs per image in (Ford, 2013).(Nurminen et al., 2015) argue that this important number of GCPs is required to tie the images to the ground reference system because camera calibration information as well as adequate approximate exterior orientations are often missing.Thus, selfcalibration is usually re-performed during the absolute image orientation stage.In the state-of-the-art, there is also no automatic method that determines the image coordinates with the ground information acquired.All authors mention this task as very time consuming.Indeed, apart from measuring their coordinates, finding these points is very challenging.They have to be stable over time on possibly a long period of time and precisely identified on the aerial photographs (Micheletti et al., 2015).The planimetric and altimetric accuracies of the absolute image orientation depends, among other parameters, on the scale (spatial resolution) of the image survey.Most studies obtain on several image surveys <1 m to 2 m accuracies, depending on the scale of the acquisitions (Micheletti et al., 2015, Aguilar et al., 2013, Cardenal et al., 2006, Necsoiu et al., 2013, Véga and St-Onge, 2008, Fox and Cziferszky, 2008) or between 2 m to 10 m (Nocerino et al., 2012).Thus, absolute image orientation remains the major lock to full automation or, at least, to reduced computing times.
Very few methods propose solutions to overcome this issue.(Korpela, 2006) adopts multi-date bundle-block adjustment with recent images whose absolute orientations are known.He found that automatic tie-points matching is not sufficient: manual multitemporal tie-points had to be identified (279 tie-points for 288 images).(Cléry et al., 2014) and (Nagarajan and Schenk, 2016) propose to use linear features characterising permanent structures such as main road axis and building outlines as candidate ground reference.Manual digitization of these line features is still mandatory both in the archival images and the ground reference information in (Nagarajan and Schenk, 2016).(Cléry et al., 2014) opted for an automatic procedure that matches objects of a topographic vector database to lines extracted within the images.However, no bundle block adjustment was performed: the method has to be applied to each archival image independently.Results are promising even though residuals do not meet accuracy requirements.

DSM and Orthoimages
As long as the image model (absolute orientation and calibration) is known, orthoimages and DSMs can be calculated.The DSM is later produced thanks to a dense matching algorithm.Thus, the DSM can then be used to orthorectify each archival image of the aerial survey: true orthoimages can then be created (Figure 2).At that step of the process, a radiometric equalisation of the images can be performed.Finally, the mosaicking step aims at merging each orthorectified image into a single orthoimage covering the whole area of interest.Very few studies tackled this last step.A DSM on archival aerial images was produced in (Nurminen et al., 2015, Cardenal et al., 2006, Nocerino et al., 2012, Véga and St-Onge, 2008, Walstra et al., 2004, Giordano et al., 2017).The important number of GCPs used in these studies seem to allow to obtain consistent altimetric DSM accuracies: from a metric accuracy (Nurminen et al., 2015, Cardenal et al., 2006) to 4 m (Walstra et al., 2004).This altimetric accuracy depends among other parameters on the scale of the aerial photographs.However, when the image model is sufficiently good, dense image matching algorithm can be easily used on archival stereoscopic aerial image surveys.
Only (Nebiker et al., 2014) compare 3 methods for dense image matching on different criteria.In contrast, issues that are related to the generation of orthoimages (equalisation of the radiometry of a set of images, vignetting correction (Kim et al., 2010)) are rarely considered.Equalisation is often embedded in dense matching solutions: radiometric differences are minimised during the mosaicking process.However, very few models have been implemented for archival images.

Proposed approach and contributions
Series of archival aerial images have no counterpart for the computation of long-term and relatively dense sampling of the Earth surface.This data can chronicle planimetric and altimetric changes that could have never yet been studied.However, in the current state of knowledge, it appears that the determination of sufficient ground reference is the major obstacle to process larger areas and more epochs.This statement concerns both non-photogrammetric and photogrammetric approaches.It can as well be noticed that very few studies propose to automate this stage.Therefore, we propose an algorithm that brings solutions to address the problems mentioned above.Our method presents several contributions: 1.The automatic determination of dense ground references, stable over time, in multi-view images; 2. A coarse-to-fine strategy that benefits from the pre-computed orthoimages and Digital Surface Model so as to refine the orientation on the historical surveys on a recent master dataset.
In this paper, a photogrammetric approach is adopted considering that the DSM that can be computed is an important asset to reach full automation.In the archival aerial image processing, automatic tie-point extraction and matching has never been tackled between two epochs when land-cover changes are important (time interval) or when the images of the surveys are very different (in terms of radiometry and spatial resolution).This issue is overcome with such a two-step absolute image orientation strategy.First, coarsely georeferenced orthoimages and DSM are produced for each year without any additional ground reference information (except basic metadata, usually stored with the archival aerial images).Then, a fine absolute image orientation of the archives is computed.The determination of ground references is performed using the standard tie-point extraction + matching solution, between the coarse orthoimages and a recent orthoimage/DSM (GCPs ground coordinates measurement).The coarse absolute image orientation and coarse DSM are inserted so as to automate the identification of the selected ground references in the archival aerial images (GCP image coordinates).
In practice, MicMac software (Pierrot Deseilligny and Cléry, 2011) for used for all steps except except automatic GCP extraction.

METHOD
Figure 3 shows the photogrammetric workflow proposed in this paper with main inputs and outputs.The workflow can be divided into 3 main parts.First, the interior and relative orientation steps are based on a standard state-of-the-art photogrammetric method (Section 3.2).Then, an original two-step (coarse-to-fine) absolute image orientation is proposed (Sections 3.3 and 3.5).Eventually, the automatic GCPs determination phase provides ground references for a new bundle adjustment procedure (Section 3.4).

Flight plan
Orthoimage DSM

Interior and relative orientation
Interior orientation

Inputs
The necessary inputs of the method are (i) the scanned archival photographs (ii) accompanied with very basic metadata information.The metadata can correspond to the digitization process (scan resolution) or the configuration of the aerial survey (focal length of the camera, the physical size (mm) of the photographs and a coarse flight plan).Focal length and the physical size of the photographs are metadata that had usually been archived.A coarse flight plan is also necessary, albeit it does not need to be very accurate.The planimetric ground coordinates of the image nadir (≈ 100 m) are already available in most cases.They are bound to be computed since a basic requirement for the distribution of the images though web-based platforms.Approximate flight height is also necessary in order to have a coarse 3D information on ground position of the pose of the images.This information is commonly written on the photographs.Recent orthoimages and DSMs are also required as a planimetric and altimetric reference.They act as a master dataset on which slave archival images are oriented, date-by-date.

Interior and relative orientation
For interior and relative orientations, state-of-the-art methods usually implemented in photogrammetric softwares are used.The literature shows that these existing methods are suitable to process each archival aerial survey independently.First, the interior orientation is performed with a semi-automatic method implemented in MicMac software (Pierrot Deseilligny and Cléry, 2011).The image coordinates of the fiducial marks are determined on a single archival image and are thereafter found automatically on the other images of the aerial survey using pattern matching techniques.Archival images are resampled after interior orientation.A common relative orientation and a camera selfcalibration are performed through a least square bundle-block adjustment based on automatic tie-point extraction and matching.At this stage a first camera calibration is thus estimated.

Coarse absolute orientation
The coarse absolute orientation proposed in this paper consists in applying a geometric transform (3D similarity) to the relative image orientations.This transform is estimated comparing the archival image centre ground coordinates of the relative poses and the flight plan.Coarse absolute image orientation is then obtained.Next, orthoimages and DSM are computed using this coarse image orientation.

Automatic ground reference identification
Several approaches have already been proposed in order to perform cross-domain homologous points extraction out of multimodal remote sensing images as for instance using Dense Adaptive Self-Correlation descriptors (Kim et al., 2015).Methods have also been developed so as to deal with diachronic images.They have never been applied to archival aerial images.They are often learning approaches: keypoint detection and matching are then trained from reference examples.Among these methods, TILDE (Verdié et al., 2015) formulates this problem as a regression model while LIFT (Yi et al., 2016) adopts deep learning models.The approach proposed in (Aubry et al., 2014) is adopted here.Indeed, it has been specifically developed in order to extract tiepoints between images and historical paintings.Thus, it can deal with both multi-modal and multi-temporal data.This method relies on HoG descriptors (Dalal and Triggs, 2005) and also considers matching as a classification problem.The underlying idea is to train for each pixel's HoG descriptor a binary classifier "similar descriptor"/ "non similar descriptor" from one positive example and many negative examples.A multi-resolution approach is adopted for keypoint detection.A scale-space is computed from the image, and keypoints are detected at each level as local maxima of the keypoint relevance score, assessing how much this point is discriminative according to its HoG descriptor.This keypoint relevance score is related to the optimal cost a squareloss SVM.In practice, keypoints are first extracted out of the recent reference orthoimages.Thus, planimetric ground coordinates are directly associated with these keypoints.Height information from the recent reference DSM can then be linked to them.This enables to retrieve 3D GCPs.They are then matched to the archival orthoimage calculated during the previous step.A tile-based approach is adopted.Indeed, the coarse georeferencing of archival data makes it possible for each tile to select only the corresponding subset of possible keypoints detected from the reference orthoimage.A dense matching is then performed: a scale-space of the archival image tile is calculated, and the matching scores between HoG descriptors are calculated (i) for each detected keypoint and (ii) for each pixel of each level of this scale-space.This matching score is defined from a Linear Discriminant Analysis classifier.For each reference keypoint, the two best matching scores are considered.The ratio between the scores of the first and second best matches is computed (Lowe, 1999).At the end, only a subset of matches with the best matching scores ratio is kept.A RANSAC filtering is then applied to them.For the remaining matches, points detected in the archival orthoimage are reprojected in the corresponding archival aerial images, to retrieve their image coordinates.This re-projection is performed using the ground-to-image model given by the coarse absolute image orientation and the height information from the archival DSM generated at the previous step.

Fine absolute orientation
Absolute archival image orientations as well as camera calibration are re-estimated at this step: the automatic GCPs determined in Section 3.4 are used together with tie-points within a robust least square bundle block adjustment.The coarse absolute image orientation is used as the initial solution.The fine absolute orthoimages and DSM are finally computed.

Experimental set-up
The method was tested on two study areas.The Fabas study area is an agricultural and forested areas where limited urban changes have occurred.The site is located in the South Western France.Important land-cover changes correspond to the regrouping of agricultural parcels, the creation of an artificial lake and natural forest evolution.The Frejus study area is located in the South-East of France.On this area, major urban land-cover changes had occurred.The site was mainly dedicated to agriculture with small villages after the second World War.It has gradually become dominated by continuous dense and residential urban areas.5 different aerial surveys for each study area were selected (Table 1).Only the panchromatic channel is kept.
In the proposed method, reference orthoimages and DSM are needed (master).They are extracted from the reference French national image databases.Orthoimages were produced at a 0.5 m GSD from airborne images.An orthoimage of 2014 was used for Frejus study site, and a one of 2010 for Fabas.DSMs were computed at a spatial resolution of 0.5 m.   4c,d shows that the scale of the coarse absolute products is consistent but the planimetric accuracy appears biased (only on the y-axis (∼120 m) on this area).Such shift can also be observed on Figure 6 for 1942 survey.The transform was found rigid as a simple first order polynomial transform applied to the coarse orthoimage with few manual GCPs found in the reference orthoimage (Figure 4d) is sufficient to coarsely correct it.The resulting corrected orthoimage is given in background orthoimage of Figure 4c).However, a comparison between the coarse absolute and the reference DSM show that the altimetric residuals are still significant and systematic.No polynomial transform was found to be able to correct these residuals.This problem is probably due to the difficulty to estimate the camera calibration only with the tie points.agery of 1954.It can be seen all extracted points do not have the same precision but even the worst of them could generally at least provide anchors to improve a coarse georeferencing.It can also be noted that GCPs were generally detected near man-made structures (roads and buildings), which could lead to inhomogeneous density in natural areas.These GCPs were then used within a robust least square bundle adjustment, allowing to filter some errors.The number of remaining GCPs as well as the mean absolute residuals are provided in Table 2. Orthoimages calculated from archival images have been first visually assessed on both areas : they were displayed together with current orthoimages in a chessboard view to assess the continuity and the alignment of persistent landscape features.They were also superimposed with recent reference building and road vector databases.Representative examples of such views are shown on Figures 7 and 8 for both study areas.The 2D registration was generally satisfactory, corresponding to what could be expected according to the quality of the different input data.At first glance, the quality of DSMs seems suitable : the different ground/above ground objects can be correctly distinguished.In addition, when considering the difference between the reference DSM and the calculated archival DSM, errors due to systematisms were strongly reduced (<2 m in unchanged areas).The difference is shown in Figure 9 for the 1971 Fabas study area.However, some problems still occur.For instance, some linear discontinuity jumps can still be observed, revealing some remaining calibration errors.Stronger gaps can be observed near the borders of the study site, where the bundle adjustment was less constrained.Nevertheless, it must be kept in mind that both study areas are not simple cases.For instance, an important part of Frejus area is covered by the sea where not tie-points and GCPs can be found, thus limiting the possible constraints during the calibration process.

Orthoimages
Table 2. Number of automatic GCPs extracted (nb) and comparison of their 3D ground coordinates (after fine absolute georeferencing vs. current reference).

CONCLUSION AND PERSPECTIVES
An almost fully automatic workflow for the photogrammetric georeferencing of images from archival aerial photogrammetric surveys was proposed based on a coarse-to-fine absolute image orientation approach.Homologous point extraction between each aerial survey and a current master reference was performed automatically using GCPs.The computation of the coarse absolute orthoimages makes this step much easier.Besides, advantage was taken from the coarse absolute image orientation and the computed DSM to perform this matching step on true orthoimage and then to automatically determine the position of the homologous point in the archival images.Archival orthoimages and DSMs could then be calculated.Their planimetric registration appeared to be very satisfactory.3D information is improved but still remains imperfect.The quantitative analysis provided in this paper is not sufficient.Indeed, it must be kept in mind that the extracted GCPs have varying geometric precision.Indeed, even for correct matches, their planimetric depends on the level of the scale-space at which they are detected.Besides, their height is estimated from a DSM.Its precision is related to the one of the DSM.Thus, some quite strong residuals can in fact correspond to the truth (land-cover changes) and not to an erroneous bundle adjustment result.A rigorous planimetric and altimetric accuracy assessment will be proposed in the near future.3D check points are intended to be measured by human operators soon.Additional improvements are still necessary, especially to more numerous and robust GCPs.First, more GCPs could be extracted using other point-based methods, as for instance deep learning based ones (Yi et al., 2016).It could be possible to train a specific architecture, benefiting from available georeferenced orthoimages at different dates.Other approaches could also be investigated.Some buildings are stable features across time.Thus, as the DSM was improved, it can now be considered to perform 3D registration between archival DSM and reference DSM over patches around them, taking into account a quality score to assess whether a change occurred.The approach of (Cléry et al., 2014) could also be used to match segments detected from the archival orthoimage to the roads and buildings from topographic databases.Such ground control lines could then be integrated in the bundle adjustment.Eventually, the planimetric accuracy obtained can be consistent enough to train land-cover classification models using recent land-cover databases.The tie-points and automatic GCPs could then be labelled with a rather coarse landcover class (e.g., vegetation/non-vegetation), and consequently filtered out for the final adjustment.

Figure 1 .
Figure 1.Spatio-temporal overview of state-of-the-art studies.Papers are plotted according to the area of the study site, the number of dates and the time range of the time series.

Figure 2 .
Figure 2. Standard photogrammetric pipeline for archival image processing.

Figure 3 .
Figure 3. Proposed method.See text for more details.
nb: number of images of the aerial survey; scale: mapping scale of the aerial survey; scan: the scan resolution during digitization; GSD: an estimation of the ground sample distance; focal: the focal length of the camera; height: the flight height of the aerial survey; size: the physical size of the archival photographs.meaningful: roads, forests and hedges are clearly visible on this natural area.The different archival images were correctly merged together to produce the coarse absolute orthoimage.No local displacements were observed between each orthorectified archival image.Comparing Figure 4a,b and Figure

Table 1 .
and DSM assessment 4.2.1 Coarse absolute orthoimages and DSMs An example of the coarse absolute orthoimages and DSMs is given in Figure 4 for a small part of the Fabas study area(1984 aerial survey).Visually, the level of details contained in the DSM is highly Aerial survey datasets.