Lucien Wald

Thierry Ranchin

Marc Mangolini

FUSION OF SATELLITE IMAGES OF DIFFERENT SPATIAL RESOLUTIONS: ASSESSING THE QUALITY OF RESULTING IMAGES (1)

Keywords: fusion, image quality, multisensor, multiresolution. Running satellite images of different spatial resolutions

Methods have been proposed to produce multispectral images with enhanced spatial resolution using one or more images of the same scene of better spatial resolution. Assuming that the main concern of the user is the quality of the transformation of the multispectral content when increasing the spatial resolution, this paper defines the properties of such enhanced multispectral images. It then proposes both a formal approach and some criteria to provide a quantitative assessment of the spectral quality of these products.

Five sets of criteria are defined. They measure the performance of a method to synthesize the radiometry in a single spectral band as well as the multispectral information when increasing the spatial resolution.

The influence of the type of landscape present in the scene upon the assessment of the quality is underlined, as well as its dependence with scale. The whole approach is illustrated by the case of a SPOT image and three different standard methods to enhance the spatial resolution.

INTRODUCTION

Many works have recognized the benefit of merging high spectral resolution (or spectral diversity) and high spatial resolution images, particularly in land mapping applications. Some satellite systems offer both kinds of data. The Landsat Thematic Mapper (TM) system presents six spectral visible and nearinfrared bands with a 30 m resolution and a thermal infrared band with a 120 m resolution. The present SPOT system has three spectral bands with a 20 m resolution (called XS1, XS2, XS3) together with a panchromatic band of 10 m resolution (called P). The next SPOT systems will have more spectral bands and enhanced spatial resolutions, and will also retain this dual resolution capability. Besides these widely used systems, other systems offer either spectral capabilities or high spatial resolution, including the Russian KFA and KVR systems, with resolutions of 10 to 2 m, and other systems to come from United States and elsewhere. This list does not pretend to be exhaustive and is intentionally restricted to the systems presently used for fine-scale land mapping. It illustrates the diversity of available data which can be merged together for a better knowledge of our environment.

Many methods have been proposed for the merging of high spectral and high spatial resolution data in order to produce multispectral images having the highest spatial resolution available within the data set.

We are only concerned with those methods which claim to provide a synthetic image close to reality when enhancing the spatial resolution, and not those which only provide a better visual representation of the image (e.g., [START_REF] Carper | The use of intensity-hue-saturation transformations for merging SPOT Panchromatic and multispectral image data[END_REF]. Some methods are specific to the case of SPOT [START_REF] Anonymous | Guide des utilisateurs de données SPOT, 3 volumes[END_REF][START_REF] Pradines | Improving SPOT images size and multispectral resolution[END_REF], while others are more general (Blanc et al., 1996;[START_REF] Chavez | Comparison of three different methods to merge multiresolution and multispectral data: Landsat TM and SPOT Panchromatic[END_REF][START_REF] Garguet-Duport | Analyse spatiale d'une zone alluviale par une nouvelle méthode de fusion d'images SPOT multispectrales (XS) et SPOT panchromatique (P)[END_REF][START_REF] Garguet-Duport | The use of multiresolution analysis and wavelets transform for merging SPOT panchromatic and multispectral image data[END_REF][START_REF] Li | Multisensor image fusion using the wavelet transform[END_REF][START_REF] Iverson | Adaptive image sharpening using multiresolution representations[END_REF][START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF][START_REF] Mangolini | Fusion d'images SPOT multispectrales (XS) et panchromatique (P), et d'images radar[END_REF][START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF][START_REF] Pellemans | Merging multispectral and panchromatic SPOT images with respect to the radiometric properties of the sensor[END_REF][START_REF] Tom | System for and method of enhancing images using multiband information[END_REF]Yocki, 1996).

These methods make use of the data having the best spatial resolution to simulate multispectral images at this resolution. However no methods propose any assessment of the quality of the resulting synthetic images except for Blanc et al., Mangolini et al., Munechika et al. and Li et al., the last approach being rather inapplicable to Earth images.

The goal of the present paper is to propose both an approach and some criteria for a quantitative assessment of the quality. In doing this, we assume that the main demand of the user concerns the quality of the transformation of the multispectral content when increasing the spatial resolution. The better the simulation of the spectral content at the enhanced resolution, the more accurate the classification for mapping purposes. Our approach is similar in essence to the one used by [START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF] and, for the case of SPOT, by [START_REF] Mangolini | Fusion d'images SPOT multispectrales (XS) et panchromatique (P), et d'images radar[END_REF][START_REF] Mangolini | Evaluation de la qualité des images multispectrales à haute résolution spatiale dérivées de SPOT. In Compte-Rendus du colloque "Qualité de l'interprétation des images de télédétection pour la cartographie[END_REF] and [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF]. The present work offers a more extensive discussion of the approach and proposes a formalism as well as quality criteria. The influence of the resolution itself with respect to the type of landscape upon the assessment of the quality is finally discussed.

BACKGROUND

At first some general comments are needed, in particular regarding the superimposability of images and the simultaneity of their acquisition by the sensors. The above-mentioned methods require superimposable images, once all images are set to the lowest available spatial resolution (e.g., 20 m in the case of SPOT). Some images are already co-registered, such as Landsat images. Otherwise this can be done by means of standard methods available in public or commercial software packages for image processing. The images of lowest resolution (e.g., XS in the case of SPOT) are projected into the geometry of the image of highest resolution (e.g., P in the SPOT case) degraded to the lowest resolution (e.g., 20 m in the SPOT case). During the process, a resampling of the multispectral images is made. The resampling operator has an influence upon the final result. In most cases, a bi-cubic interpolator offers a good compromise between the accuracy of the result and the required computer time. In the following, for the sake of the simplicity, the term "image of lowest resolution" will denote the projected resampled image of lowest resolution, if this is required.

The images of different resolutions may not have been acquired simultaneously. Optical properties of the atmosphere are different from one date to the other. This induces further spectral distortion in the set of images. As for the landscape, as long as the time-lag is small with respect to the time scale of the variations in small-size features, its influence upon the quality of the transformation of the spectral content when enhancing the spatial resolution is low or negligible. Such a time scale is greatly variable and does depend upon the objects themselves as well upon the spatial and spectral resolutions with which they are observed. If the time-lag is large, the user must weight its consequences. He should know precisely the merging method to be used, because all methods do not take into account in the same way the small structures to be injected from the images of highest resolution into the images of lowest resolutions. For example, [START_REF] Anonymous | Guide des utilisateurs de données SPOT, 3 volumes[END_REF] recommends that its method should be used only for coincident SPOT XS and P data. This method from the CNES (French space agency) is called P+XS and is certainly the most commonly used method among those cited above. However, it is frequently applied to non-coincident SPOT data because of the difficulties in obtaining coincident couples of images.

The quality of the resulting synthetic images is usually assessed by visual inspection, a necessary step. For example, Figures 1 to 4 present a SPOT sub-scene of Barcelona., acquired on 11 September 1990.

Barcelona is a large city located in north-east Spain, on the Mediterranean coast. Its harbor is the busiest in Spain. This extract of the SPOT scene displays the western, newest districts of the city. Displayed is only a part (360 by 360 pixels) of the sub-scene (512 by 512 pixels) which is dealt with in the following sections. The sub-scene is mostly comprised of urban districts, highways and railroads. It also exhibits small agricultural lots and mountainous areas covered by typical Mediterranean vegetation. Such an urban area has been selected for illustration because it is certainly the most difficult type of landscape to deal with according to our knowledge. Urban areas often point out the qualities and drawbacks of algorithms because of the high variability of information in space and spectral band, induced by the diversity of features both in size and nature.

The panchromatic band P is shown in Figure 1, while the XS1 is in Figure 2. The latter has been magnified by a factor of two. In Figure 3 is the synthetic image obtained by the CNES method (XP1), and in Figure 4, the synthetic one obtained by the ARSIS method (XS1-HR). These methods are discussed later. When comparing images, one must pay attention to the contrast table (look-up table) because it acts as a filter (together possibly with the printer) between the information and the human observer. In the case of SPOT, the radiances observed in the P, XS1, and XS2 bands are very similar for a spectrally neutral target. For a spectral band i, the radiance R i is linked to the digital count DC i by the calibration factor a i : i.e., R i = DC i / a i In this particular case, the calibration factors are very similar for the P and XS1 bands (see Table 1), and, thus, so are the digital counts. It follows that the same look-up table can be applied to each image in Figures 1 to 4 and that they can be visually compared. Beyond demonstrating the interest of merging P and XS data, visual inspection clearly shows the major drawbacks of both methods. In Figure 3, local contrasts are too much reinforced. The extreme values are also reinforced: the white areas are whiter, compared to Figure 2, and the dark areas are darker. In Figure 4, on the contrary, one may think that local contrasts are too smooth, but gray tones are very similar to Figure 2.

The objective comparison of the visual quality of multiple images is a difficult and lengthy task to handle. The human visual system is not equally sensitive to various types of distortion in an image. The perceived image quality is strongly dependent upon the observer and also upon the thematic application. Standard protocols have been defined, namely in the field of television and image compression. Such a protocol is described, for example, in [START_REF] Lu | A comparative study of wavelet image coders[END_REF] for the compression of still images. A panel of human viewers judge some well-defined aspects of the images. Then their notations are weighted and further processed to obtain a mean opinion score defining the quality of the result. When it comes to the assessment of the quality of a set of multispectral images, the mass of data becomes very large. This dramatically increases the difficulty in computing a quantitative picture quality scale. Beyond the visual inspection, mathematical criteria are needed.

One simple way is to look at the histograms of the synthetic products and to compare them to the original ones. The histograms for the previous images (Figures 1 to 4) are presented in Figure 5. On the upper half are those for the original images P and XS1. For the latter, the resolution is 20 m only: it contains four times fewer pixels than P or the synthetic images. For the comparison, the XS1 histogram has been normalized to the others by multiplying the number of pixels by four. Though the resolution is increased by a factor of two relative to XS1, the histograms of the synthetic images are expected to be close to the XS1 one in shape. This is true for the XS1-HR histogram. Its highest frequency is close to four times the peak of XS1. On the contrary, the XP1 histogram is much closer to the P one, both in shape and in peak. This comparison of histograms is a fairly good estimator of image quality, and is easy to handle. However, the effect of the spatial resolution upon the statistical properties of an image should not be neglected (e.g., [START_REF] Kong | Effet de la résolution spatiale sur les propriétés statistiques des images satellites : une étude de cas[END_REF][START_REF] Woodcock | The factor of scale in remote sensing[END_REF][START_REF] Raffy | Remotely-sensed quantification of covered areas and spatial resolution[END_REF][START_REF] Welch | The effect of spatial resolution upon the texture-based cloud field classifications[END_REF]. That means that we should not try to identify the statistical properties of a synthetic product to those of the original image. These discrepancies in statistics depend upon the observed type of landscape.

Therefore, other mathematically-sound criteria are needed.

Another approach is to compare land-use maps obtained after spectral (and possibly textural) classification of the synthetic images. These maps are compared either to the map obtained from original low-resolution data (e.g., SPOT XS), or to ground truth. In the first case, the same assumption as above is made, that some statistical properties are preserved through the increase in resolution. Hence, this approach should be avoided. More generally, classification greatly reduces the content in information, making discrepancies between methods disappear. For example, [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF] compared their method to Price's one (1987) for Landsat-TM. On the one hand they obtained a large relative difference between both methods -which we have estimated at about 25 % (RMS) according to the figures available in their paper -while, on the other hand, the classification rates were very similar. This classification approach is valuable because land-use mapping is often the goal of satellite image processing. However, it may not reflect the overall performance of a method since the results depend too much upon the type of landscape, its diversity, its heterogeneity, the time of observation, the optical properties of the atmosphere, the sensor system itself (including the viewing geometry), the type of classifier (supervised, unsupervised), and the classifier itself.

The type of landscape present within the image used to assess the quality of a synthesizing method has a strong influence upon the results. Whatever the method, the more predictable the change in signal with the scale, the better the quality of the final product. Over areas such as the ocean or large agricultural lots, which appear very homogeneous at, say, 20 m resolution, the error made in assuming that these areas are still homogeneous at, say, 10 m resolution, is small. On the other hand, urban areas or small agricultural lots are among the most difficult cases because they exhibit a large number of interwoven objects having different scales. [START_REF] Ranchin | Fusion of images and raster-maps of different spatial resolutions by encrustation: an improved approach[END_REF] examined the SPOT image of Barcelona presented above.

They found that, for the homogeneous part covering the Mediterranean Sea (not visible in Figures 1 to 4, and actually not within the sub-scene used for illustration), all information, expressed as variance, is borne by structures larger than 40 m. On the contrary, for the urban area, half the information is borne by structures having sizes less than 40 m. Such cases do not possess self-similarity properties, though some parameters, such as the growth of city limits, can be approximated by fractal functions (e.g., [START_REF] Batty | Cities as fractals: simulating growth and form[END_REF].

In other words, structures observed at, say, 10 m resolution, cannot be accurately predicted from their observations at lower resolution, say, 20 m. This is well-known by experienced image interpreters, and is also sustained by mathematical evidence (e.g., [START_REF] Fung | Spatial composition of spectral classes: a structural approach for image analysis of heterogeneous land-use and land-cover type[END_REF][START_REF] Ramstein | Structures spatiales irrégulières dans les images de télédétection[END_REF][START_REF] Woodcock | The factor of scale in remote sensing[END_REF]. The benefit of an image of a higher spatial resolution is the greatest in these cases. Hence, we recommend that test images should mainly include such areas. Such cases also offer a large diversity of spectral signatures, which is helpful in judging the ability of a method to synthesize the spectral signatures during the change in spatial resolution.

A FORMAL APPROACH FOR QUALITY ASSESSMENT OF THE RESULTING SYNTHETIC

IMAGES

Let denote the images of lowest resolution by M l , and the images of highest resolution by P h . The subscripts l and h denote the resolution of image M or P, i.e., low and high resolution, respectively. The subscript s (for small) is to be used later: it denotes a resolution which is lower than l, say, e.g. two times lower. According to the previous section, the images M l and P l are superimposable. Both M l and P h have been obtained by a sensor. The merging method aims at constructing synthetic images M* h .

These synthetic images must have the three following properties. First, any synthetic image M* h , once degraded to its original resolution l, should be as identical as possible to the original image M l . For example, in the case of SPOT, the synthetic image is called XS* 10 . Once resampled to 20 m, this image should be as close as possible to the original XS image. Besides the effects induced by time-lag in image acquisition, as discussed in previous section, observing this property means that the method takes into account the differences in atmospheric effects affecting the images of lowest and highest resolutions which have not been acquired within the same spectral bands. Second, any synthetic image M* h should be as identical as possible to the image M h that the corresponding sensor would observe with the highest resolution h. Because the synthetic image M* h do not match M h exactly, a property should be added, dealing with the entire set of channels. Third, the multispectral set of synthetic images M* h should be as identical as possible to the multispectral set of images M h that the corresponding sensor would observe with the highest resolution h. In the case of SPOT, this set of synthetic images is the triplet (XS1* 10 , XS2* 10 , XS3* 10 ). The assessment of the quality of the resulting spatially enhanced spectral images is now equivalent to the verification of these properties.

Testing the first property: Any synthetic image M* h , once degraded to its original resolution l, should be as identical as possible to the original image M l . To achieve this, the synthetic image M* h is spatially degraded to an approximate solution M' l of M l . If the first property is true, then M' l is very close to M l .

The difference between both images is computed on a per-pixel basis. This difference image should be visually compared to the original image in order to detect trends of error, if any, possibly related to the type of landscape. Then some statistical quantities are to be used to quantitatively express the discrepancies between both images. These quantities are similar to the first and second sets of criteria described under the second property below.

An important point here is the way the synthetic image M* h is degraded to M' l . Some wavelet transforms have the ability to separate scales well, that is, to separate small size structures from larger ones and, therefore, to fairly well simulate what would be observed by a lower resolution sensor (e.g., Ranchin and Wald, 1993). [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF] used an averaging operator on a window of 3 by 3 pixels. Such an operator does not have this ability and is not as appropriate here. Other filtering operators can be used, some of them simulating a given modulation transfer function (MTF) of a sensor. A comparison was made on a few scenes using some operators, such as the truncated Shannon function, bi-cubic spline, pyramid-shaped weighted average, and wavelet transforms of Daubechies (1988, regularity of 2, 10 and 20). It showed relative discrepancies between the results on the order of a very few per cent. In conclusion, there is an influence of the filtering operator upon the results, but it can be kept very small provided the operator is appropriate enough.

Testing the second property: Any synthetic image M* h should be as identical as possible to the image M h that the corresponding sensor would observe with the highest resolution h. The second and third properties are difficult to test, because they refer to M h , an image that would be sensed if the sensor had a better resolution. This image, of course, is not available; otherwise, all the above-cited methods would not have been developed. The difficulty is partly overcome by the following approach:

-The available images P h and M l are degraded to respectively P l and M s , respectively. For the SPOT case, the P and XS images are degraded to P 20 (20 m resolution) and XS 40 (40 m resolution). The images P l and M s are very close to what the corresponding sensor would have measured with a degraded resolution, as discussed previously.

-Then the synthesizing method under assessment is applied to P l and M s . It provides a synthetic image M* l (XS* 20 in the SPOT case).

-This synthetic image M* l is compared to the image-truth M l (XS in the SPOT case) by means of some criteria described below. The numerical comparison should be made preferably in physical units and also in relative values. Thus, different tests made over different scenes may be compared.

-This comparison provides an assessment of the quality of M* l . It is assumed that this quality is fairly similar to that of the synthesized high-resolution image M* h . This point will be discussed later.

Such an approach was made by [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF] and [START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF][START_REF] Mangolini | Evaluation de la qualité des images multispectrales à haute résolution spatiale dérivées de SPOT. In Compte-Rendus du colloque "Qualité de l'interprétation des images de télédétection pour la cartographie[END_REF]. To assess the quality of M* l , a difference is computed between M l and M* l in a way similar to that used for the first property. After visual inspection, the difference image is reduced to a few statistical parameters which summarize it. There are a large number of candidate parameters. We have computed many for several tens of cases. We have retained some whose definitions are well-known to engineers and researchers and which clearly characterize the advantages and disadvantages of a method.

Two sets of criteria are proposed to quantitatively summarize the performance of a method in synthesizing an image in one spectral band. The first set of criteria provides a global view of the discrepancies between the original image M l and the synthetic one M* l . It contains:

-The bias, as well as its value relative to the mean value of the original image. Recall that the bias is the difference between the means of the original image and of the synthetic image. Ideally, the bias should be null.

-The difference in variances (variance of the original image minus variance of the synthetic image), as well as its value relative to the variance of the original image. This difference expresses the quantity of information added or lost during the enhancement of the spatial resolution. For a method providing too many innovations (in the sense of information theory), i.e., "inventing" too much information, the difference will be negative because the variance of the synthetic image will be larger than the original variance. In the opposite case, the difference will be positive. In information theory, the entropy describes the quantity of information. However, we selected the variance difference because most researchers and engineers are much more familiar with variance, and entropy and variance act quite similarly for our purpose. Ideally, the variance difference should be null.

-The correlation coefficient between the original and synthetic images. It shows the similarity in small size structures between the original and synthetic images. It should be as close as possible to 1.

-The standard deviation of the difference image, as well as its value relative to the mean of the original image. It globally indicates the level of error at any pixel. Ideally, it should be null.

The error at pixel level may be more detailed. Let us compute at each pixel the absolute relative error (the absolute value of the difference between the original and synthetic values, divided by the original value).

Then the histogram of these relative errors is computed. It can be seen as the probability density function.

Therefore, we can compute the probability of having at a pixel a relative error (in absolute value) lower than a given threshold. This probability denotes the error made at pixel level, and hence indicates the capability of a method to synthesize the small size structures. The closer to 100 percent the probability for a given error threshold, the better the synthesis. The ideal value is a probability of 100 percent for a null relative error. Here, for reasons of computer precision, the lowest threshold "no relative error or null error" is set to 0.001 percent Testing the third property: The multispectral set of synthetic images M* h should be as identical as possible to the multispectral set of images M h that the corresponding sensor would observe with the highest resolution h. Visual inspection may be made through color composites of, for example, the first three principal components of the set of images. Both color composites should agree visually. Most methods for color composites are using dynamical adjustment for color coding (e.g., [START_REF] Albuisson | Codage trichrome et classification[END_REF]. If the sets of images are different, even slightly, then the color coding will be different for both composites and no comparison will be possible. Practically, we recommend the following approach. For each spectral band, the M l and M* l images are juxtaposed into a single computer file. The principal components analysis as well as the color coding are performed on this set of files. The projected M l and M* l images are then extracted from these projected files and the color composites are displayed, simultaneously or alternatively, onto the screen. This approach guarantees that the color composites are comparable. Of course, if only three spectral bands are available as in the SPOT case, there is no need to perform a principal components analysis. The advantage of this visual assessment is that it does show trend in errors, if any, possibly related to landscape features. The drawback of it is that it is a subjective assessment and also that this assessment may be limited either by physiological factors (e.g., color contrast perception by humans), or by technical factors (e.g., when a large number of spectral bands are present). In the latter case, and if the landscape offers a large variety of objects, the color re-coding of the first three principal components reduces dramatically the differences between the M l and M* l images, particularly if these differences are random, i.e., not related to a peculiar landscape feature or to a spectral band.

A quantitative assessment can be made using the following three additional sets of criteria which quantify the performance of a method to synthesize the spectral signatures during the change in spatial resolution.

The third set (numbered after the two sets described above for the second property) deals with the information correlation between the different spectral images taken two at a time. This dependence can be expressed by the correlation coefficients, with the ideal values being given by the set of original images M l . As an example, for the case of SPOT, the correlation coefficient between P 20 and XS1* 20 is computed and compared to the correlation coefficient for P 20 and XS1 20 . This is done for every pair. The fourth set of criteria partly quantifies the synthesis of the actual multispectral n-tuplets by a method, where n-tuplet means the vector composed by each of the n spectral bands at a pixel. It comprises the number of different n-tuplets (i.e., the number of spectra) observed in the original M l and in the synthesized M* l sets of images, as well as the difference between these numbers. A positive difference means that the synthesized images do not present enough n-tuplets; a negative difference means too many spectral innovations.

The previous criteria do not guarantee that the synthesized n-tuplets are the same as in the original image M l . The fifth and final set of criteria assesses the performance in synthesizing the actual n-tuplets. It deals with the most frequent n-tuplets, because they are predominant in multispectral classification. For a given threshold in frequency, only the n-tuplets having a frequency (relative number of pixels) greater than this threshold are used. The threshold is set to 0.01 percent, 0.05 percent, 0.1 percent, and 0.5 percent, successively. The greater the threshold, the lower the number of n-tuplets, but the greater the number of pixels exhibiting one of these n-tuplets. For each of the n-tuplets, the difference is computed between the original frequency and the one observed in the synthesized images. These differences are summarized by the following quantities:

-the number of actual n-tuplets, the number of coincident n-tuplets in the synthesized images, and the difference between these numbers, expressed in absolute and relative terms, -the number of pixels in these n-tuplets, in absolute and relative terms, and -the difference between the above number of pixels for the original and synthesized images, in absolute and relative terms. [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF] partly quantify the performances in synthesizing the multispectral information by first computing the root mean square (RMS) of the differences, pixel per pixel, between the synthesized M* l and original M l images, for each spectral band, and then by summing up these spectral RMS values to obtain a global error, which should be as low as possible. This global error can easily be computed from our first set of parameters, i.e., the bias and the standard deviation. Other criteria may be further defined, dealing for example with the performances in synthesizing the actual n-tuplet at a given pixel. We do not present them because, in the classification process, pixel spectral values are aggregated with their spectral neighbors. Hence, a small difference between the synthesized and the actual n-tuplets at a given pixel may have an impact on classification ranging from null to significant.

Our approach is now illustrated by an example, with emphasis on the quantitative criteria. This example consists of the same SPOT sub-scene of Barcelona as above. Obviously, our objective is the discussion of the quality assessment and not the evaluation of one or more particular methods. However, we felt it necessary to use several methods in our example for a better description of the approach. We selected a very crude method, a standard one, and an advanced one. The discrepancies between the results of these methods will be used to illustrate the advantages and limits of the criteria. [START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF][START_REF] Mangolini | Fusion d'images SPOT multispectrales (XS) et panchromatique (P), et d'images radar[END_REF]; see also Ranchin et al., 1993[START_REF] Ranchin | Efficient data fusion using wavelet transform: the case of SPOT satellite images[END_REF]Ranchin et al., , 1996a) ) developed the most advanced concept for the type of fusion of concern. They called it ARSIS, an acronym for its French name "accroissement de la résolution spatiale par injection de structures". In the ARSIS concept lies a model which permits the synthesis of missing structures of small size given the set of spectral high and low resolution data. A large variety of models can be developed based on differences in the application of the ARSIS concept by the different authors (see also Blanc et al., 1996;Ranchin and Wald, 1996;Ranchin et al., 1996b). [START_REF] Garguet-Duport | Analyse spatiale d'une zone alluviale par une nouvelle méthode de fusion d'images SPOT multispectrales (XS) et SPOT panchromatique (P)[END_REF][START_REF] Yocky | Multiresolution wavelet decomposition image merger of Landsat Thematic Mapper and SPOT panchromatic data[END_REF] and [START_REF] Yocky | Multiresolution wavelet decomposition image merger of Landsat Thematic Mapper and SPOT panchromatic data[END_REF] make use of this concept but with models of poor performance, as demonstrated by [START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF] or [START_REF] Mangolini | Apport de la fusion d'images satellitaires multicapteurs au niveau pixel en télédétection et photo-interprétation[END_REF]. The model of [START_REF] Iverson | Adaptive image sharpening using multiresolution representations[END_REF] is based upon a neural network. Here we use the most recent work made for the SPOT case by Ranchin et al. (1993[START_REF] Ranchin | Efficient data fusion using wavelet transform: the case of SPOT satellite images[END_REF], which has been selected by SPOT-Image and other French organizations.

APPLICATION TO A SPOT SUB-SCENE OF BARCELONA (SPAIN)

Two other methods are also used for illustration. One is the duplication of pixels: each original pixel, say at 20 m, provides four new pixels at 10 m resolution, each of them having the same value as their parent pixel. Though it does not take advantage of the presence of the image of higher resolution (P in the SPOT case), this method is widely used because of its simplicity. Of course, the visual aspects of such synthesized images are very bad. We use this method here as a baseline to demonstrate the benefits of more sophisticated merging methods. In the following tables, this method is denoted "dup".

The P+XS method of CNES takes into account the modulation transfer function and the spectral filter of each band P and XS. It should be applied to images acquired at the same time. However, its mathematical expression is so simple that many people are using it even for non-coincident dates. Strictly speaking, it cannot be used as such in our scheme when resolutions are degraded from 20 m to 40 m. Accordingly, we define a CNES-like method, dubbed M2 (for second method). For the synthesis of images at 10 m, the M2 method is identical to the P+XS one.

The tests described above are now applied to each method for the SPOT sub-scene of Barcelona. For the sake of clarity and to avoid duplication of tables, the example will only deal with the second and third properties. In fact, the duplication and ARSIS methods are inherently built to satisfy the first property, with reservations regarding the degradation process as discussed earlier. On the other hand, the M2 method is less satisfying. From its equations, it is obvious that there is a strong radiometric influence of P and XS2 (respectively XS1) in the synthetic image in Band 1 (respectively in Band 2). This influence does not disappear when reducing the resolution to 20 m.

In this example, the P and XS images are degraded to a resolution of 20 and 40 m respectively. Then, images are synthesized at a 20-m resolution and compared to the original XS images. Table 1 provides the mean and standard deviation of the radiances for each spectral band as well as the calibration coefficients.

Tables 2 and3 provide a global view on the discrepancies between XS and XS* 20 for each spectral band (testing the second property). The bias is null, or close, for all methods. The duplication method leads to a decrease in the quantity of information. There is no innovation brought in increasing scale. This criterion agrees with the poor visual aspects of the duplicated images. The M2 method invents too much information (up to 35 percent). The panchromatic band strongly reinforces the local contrast in the first two spectral bands. It also increases the extreme values, as already shown in the discussion of Figures 1 to 4 and of the histograms (Figure 5). The ARSIS method is closer to the ideal values but lacks innovation.

The correlation coefficient denotes the similarity between structures. It is high for all methods. It is the lowest for duplication because the latter does not invent the structures of smaller size. The standard deviation of the differences is weak. The worst result is offered by the M2 method. This criterion does not reflect the visual quality of the synthetic image, and may offer too much of a global view.

Error at the pixel level can be described better through the second set of criteria (Table 4). For the ARSIS method, almost all pixels (99 percent for XS1 and XS2, 95 percent for XS3) exhibit a relative error less than or equal to 10 percent. For the two other methods, similar results are attained at a threshold of 20 percent. Also apparent in the ARSIS method is the very high percentage of pixels exhibiting null relative errors: more than 25 percent for XS1 and XS2. The M2 method provides the worst results. This Table fully demonstrates the discrepancies that can exist between visual / qualitative and quantitative assessments of the quality of a synthetic image.

The third property deals with the multispectral character of the data set and is verified through the third to fifth sets of criteria. Table 5 shows that the M2 method increases the correlation between the spectral bands XS1, XS2 and P, as evidenced by the equations of the method. The duplication method exhibits too low coefficients, due to the lack in innovation. This table provides a first indication that the multispectral character of the synthetic images provided by the M2 and duplication methods may be only partly verified.

Table 6 reports some spectral characteristics of the scene, i.e., the total number of pixels in the image, the number of spectral triplets and the average number of pixels per triplet (ratio of total number of pixels to the number of triplets). The spectral homogeneity is defined as the inverse of the number of triplets, and is expressed in percent. It characterizes the spectral diversity of a scene: the greater this parameter, the lower the diversity. For a scene exhibiting an unique spectral object (i.e., only one triplet), this spectral homogeneity would be 100 percent. The present scene has a value of 0.002 percent, which demonstrates its spectral diversity. The performance in synthesizing the multispectral information is partly shown in Table 7, which presents the difference between the actual number of triplets and the number found in the synthesized images. The duplication method does not invent enough triplets: the number of synthesized triplets is half the actual number. The method M2 invents too many, and the ARSIS method is closer.

Actually, most of the missing or superfluous triplets have a low frequency, i.e., each of them is carried by a few pixels. This is demonstrated in Table 8, which exhibits the performances in synthesizing the most frequent actual triplets. Each triplet under consideration has a frequency of at least 0.01 percent, which corresponds to 26 pixels in this example. These 1,549 triplets have a cumulative frequency of 22 percent; that is, 57,096 pixels among 262,144 (total number) carry one of these triplets. Hence, synthesizing them is of primary importance for classification purposes. The duplication and ARSIS methods retrieve all these triplets as well as their frequencies (the difference in number of pixels, relative to the original, is low). Though the M2 method synthesizes the right triplets, it does not produce the right frequencies by far (41 percent too low). These observations are in full agreement with previous tables, and particularly Tables 2 to 4.

As said previously, the total error proposed by [START_REF] Munechika | Resolution enhancement of multispectral image data to improve classification accuracy[END_REF] can be computed from Tables 2 and3. It is found (in radiance units): i.e., duplication: 12 M2: 11 ARSIS: 7

This global measure agrees with the previous conclusions, though it might be difficult from this unique criterion to decide whether the duplication or the M2 method is more suitable on a case by case basis. For example, if a very low relative error at a pixel is requested or if preservation of multispectral content is at stake, it is obvious from Table 4 (left columns) or Table 8 that duplication should be preferred.

EXTRAPOLATION OF QUALITY ASSESSMENT TO THE HIGHEST RESOLUTION

The verification of the second and third properties of the synthetic images has been made on degraded images (e.g., in the SPOT case, we have synthesized multispectral images at a resolution of 20 m). Such an approach alleviates the lack of "truth" images. How can the assessment of quality of the synthetic images be made at the highest resolution (e.g., 10 m in the SPOT case) based upon that made at the lowest resolution (e.g., 20 m in the SPOT case) ? In other words, how can one extrapolate the quality assessment made at the lowest resolution to the highest resolution? Intuitively, one thinks that, except for objects having a size much larger than the resolution, the error should increase with the resolution, since the complexity of a scene increases with the resolution. That is, one may expect the error made at the highest resolution to be greater than that at the lowest resolution. However, several recent works have demonstrated the influence of the resolution on the quantification of parameters extracted from satellite imagery. Many works have dealt with clouds. Of particular interest are the works of [START_REF] Welch | The effect of spatial resolution upon the texture-based cloud field classifications[END_REF] for satellite imagery and of [START_REF] Kristjansson | Cloud parametrization at different horizontal resolutions[END_REF], who addresses the problem of resolution in weather prediction and climate models. [START_REF] Rowe | Incorporating landscape heterogeneity in land surface albedo models[END_REF] studies the influence of the distribution of the elementary reflectors within the pixel upon the observed signal. Also relevant are the works of [START_REF] Kong | Effet de la résolution spatiale sur les propriétés statistiques des images satellites : une étude de cas[END_REF] and [START_REF] Woodcock | The factor of scale in remote sensing[END_REF]. [START_REF] Raffy | Remotely-sensed quantification of covered areas and spatial resolution[END_REF] sets up the mathematical fundamentals to explain such a behavior in rather simple cases. The results he obtained are very similar to the ones displayed by Welch et al. (1989, Figure 3; see also Lillesand and Kiefer, 1994, Figure 7.53). All these studies demonstrate that the quality of the assessment of a parameter is an unpredictable function of the resolution.

It follows that we cannot predict the quality of the synthetic images at the highest resolution (e. g., 10 m in the SPOT case) from the assessments made with synthetic images at the lowest resolution (e. g., 20 m in the SPOT case). To illustrate this discussion, we have assessed the quality of a SPOT image synthesized at 40 m, starting from a P image degraded to 40 m and a XS image degraded to 80 m. The scene is the same as before; the method used is the ARSIS one. The results are presented in Tables 9 and10, under the heading '40 m'. In these tables are reported the results obtained for 20 m, output from Tables 2 and4.

One can see that, for all parameters, the values displayed for 20 m are better than for 40 m. Hence, the method provides better estimates in synthesizing images at 20 m resolution than at 40 m. Such comparisons were made for the three methods above-mentioned and for a few different scenes, comprising mostly urban areas. It has been found in each case that the quality was best at 20 m, and also that the ranking of a method relative to the others was the same at 20 m and 40 m. Though our conclusions were always the same, it does not prove that estimates should be better at 10 m than at 20 m.

However, we can reasonably assume that the quality of the synthetic images at the highest resolution (e.g., 10 m) is close to that at the lowest resolution (e. g., 20 m).

CONCLUSION

Many methods have been proposed for the merging of high spectral and high spatial resolution data in order to produce multispectral images having the highest spatial resolution available within the data set.

Very few propose an assessment of the quality of the resulting synthetic images. The present work proposes both a formal approach and some criteria to provide a quantitative assessment of the synthetic images. The approach is based upon simple concepts, easy to understand and easy to implement and use.

Together with the visual evaluation of the synthetic images, these criteria may be used by to select a method among others, according to its performance for the criteria which are the most important for the application. For mathematicians, these criteria, somewhat extended to more complex statistical quantities, provide a tool to assess the merits and drawbacks of a method under development.

This text is illustrated by the case of SPOT images, but the concept can be applied to other combinations of sensors. For example [START_REF] Mangolini | Procédé et dispositif pour augmenter la résolution spatiale d'images à partir d'autres images de meilleure résolution spatiale[END_REF] have assessed the performances of the ARSIS method to synthesize Landsat TM6 (thermal infrared band) at a resolution of 30 m. Table 10. Probability (in percent) for having in a pixel a relative error less than or equal to the thresholds noted in the first row. The ideal value is 100 as early as the first threshold 0.001 percent. The relative errors are in absolute value and in percent. For XS1.

  It synthesizes XS* M2 images (in radiances) at a 20 m resolution in the following way: XS1* M2 = 2 P 20 XS1 40 / [ XS1 40 + XS2 40 ] XS2* M2 = 2 P 20 XS2 40 / [ XS1 40 + XS2 40 ] XS3* M2 = duplication of XS3 40

Figure 1 .

 1 Figure 1. SPOT P image (see text for explanations). Copyright CNES -SPOT Image (1990)

Figure 3 .

 3 Figure 3. XP1 synthetic image (see text for explanations).

Figure 4 .Figure 5 .

 45 Figure 4. As Figure 3 but XS1-HR.

Table 8 .

 8 Performance in synthesizing the multispectral information. Difference between the actual frequency of a triplet (XS1, XS2, XS3) and its estimate. Only the most frequent triplets are taken into account. The total of pixels they represent amounts to 22 percent of number of pixels in the image. Each triplet has a frequency of at least 26 pixels, i.e., at least 0.01 percent of the total number of pixels.

	original	dup	M2	ARSIS

Table 9 .

 9 Some statistics on the differences between the original and synthesized images (in radiance or relative value) for XS1 band.

		0.001	1	2	5	10	20	50
	20 m	27	28	63	92	99	100	100
	40 m	23	24	56	88	99	100	100

ACKNOWLEDGMENTS

This work has been partly supported by the French Ministry of Defense (Cellule d'Etudes en Géographie Numérique, Centre Technique des Moyens d'Essais de la Délégation Générale pour l'Armement) and by CNES (Centre National d'Etudes Spatiales). The advice of Patrick Dorlet, Philippe Munier, Claude Penicand and Thierry Rousselin and of the referees is gratefully acknowledged.


Eric Guiot
email: guiot.eric_1@yahoo.fr

Charged particles trajectories under time varying magnetic fields

Keywords:

Charged particles trajectories under time varying magnetic fields

-Introduction

The motion of classical charged particles undergoing in time-varying magnetic field is of paramount interest to many areas of the modern Physics, for example in Plasma Physics. However, it seems it doesn't exist at this day analytical solutions to describes correctly these trajectories. Theories and Numerical simulations which are presented are generally based on the assumption that the magnetic moment is approximately conserved, or choosing an arbitrary choice of the expression of the potential vector.

Paper presents our contribution to this problem. As it is generally the case, the model we detail here is based on corpus of classical Physics and electrodynamics, but is obtained exploring a way which seems original. In particular paper provides coherent expressions of the potential vector and of the induced electric force, such these theoretical results can constitute significant progress in the topic.

Moreover, using these results, we present the establishment of the possible trajectories of a charged point particle under a classical time varying magnetic field, given with 𝐵 = 𝐵 0 + 𝐵 1 𝑒 𝑗𝑤𝑡 Indeed, solving the differential equations of motion we have obtained, we present a family of curves definite by combinations of Whittetaker and Lambert functions. In this part we draw certain of these original curves, and we highlight a phenomenon of resonance, between the frequency of Larmor linked to the constant part 𝐵 0 , and the pulsation 𝑤.

Note that, to simplify the problem, we have only considered the case where the varying magnetic field is applied perpendicular to the plane where the particle moves. Several possible later works are evoked in the discussion, in the field of applied and experimental Physics. Moreover, and to conclude this introduction, we should remind the reader that this work is entirely based on the laws of classical Physics and electromagnetism and that we have consequently not taken into consideration relativistic or quantum effects.

We introduce a Galilean and cylindrical system of coordinate definite with

(𝑂; 𝑒 𝑟 ⃗⃗⃗ ; 𝑒 𝜃 ⃗⃗⃗⃗ ; 𝑒 𝑧 ⃗⃗⃗ )
where 𝑂 is the origin of this system and 𝑟 the radial distance to the origin with the relation 𝐹𝑀 ⃗⃗⃗⃗⃗⃗ = 𝑟𝑒 𝑟 ⃗⃗⃗

A time varying magnetics field is applied, along the (𝑂; 𝑒 𝑧 ⃗⃗⃗ ) axis, such

𝐵(𝑡) ⃗⃗⃗⃗⃗⃗⃗⃗ = 𝐵𝑧

Under the above conditions the trajectory of the particle is planar. The induced electric field accompanying the time varying magnetic field is given by the equation of Faraday [START_REF] Jackson | Classical Electrodynamics[END_REF]

∇ × 𝐸 ⃗ = - 𝑑 𝑑𝑡 𝐵 ⃗ ( 1
)
Point particle is thus submitted to the force of Lorentz.

𝐹 = 𝑞𝐸 ⃗ + 𝑞𝑉 ⃗ × 𝐵 ⃗

To symplify the writing we decompose it as 𝐹 = 𝐹 𝐸 ⃗⃗⃗⃗ + 𝐹 𝑀 ⃗⃗⃗⃗⃗

Where 𝐹 𝑀 ⃗⃗⃗⃗⃗ is the magnetic force and 𝐹 𝐸 ⃗⃗⃗⃗ the electric. The potential vector 𝐴 is simply linked to 𝐸 ⃗ and to the electric potential 𝜑 with

𝐸 ⃗ = -∇φ - 𝛿 𝛿𝑡 𝐴 (2)
And to the magnetic field with ∇ × 𝐴 = 𝐵(𝑡) ⃗⃗⃗⃗⃗⃗⃗⃗

Moreover, as usual in this kind of problem [START_REF]Trajectory under harmonic potential and magnetic force -Archive ouverte HAL[END_REF], we assume that 𝐴 and 𝐸 ⃗ are functions of the time and the position vector 𝑟. Moreover, in our case electric force isn't the derivative of a potential, and we obtain thus the system of equations

{ 𝐴 = 𝐴(𝑟, 𝑡) ⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗ (4.1) 𝜑 = 0 (4.2) 𝐸 ⃗ = 𝐸(𝑟, 𝑡) ⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗⃗ (4.3) }

-Law of motion

We know that in our system of coordinate the speeds and accelerations are given with

{ 𝑉 ⃗ = 𝑟ė 𝑟 ⃗⃗⃗ + 𝑟𝜃 ̇𝑒𝜃 ⃗⃗⃗⃗ 𝑎 = (𝑟̈-𝑟𝜃 ̇2)𝑒 𝑟 ⃗⃗⃗ + (𝑟𝜃 ̈+ 2𝑟θ ̇)𝑒 𝜃 ⃗⃗⃗⃗ }
Which leads to obtain finally, using the principle of the dynamics a system of equations

{ 𝑚(𝑟̈-𝑟𝜃 ̇2) = 𝑞𝐵𝑟𝜃 ̇+ 𝐹 𝐸𝑟 (5.1) 𝑚(𝑟𝜃 ̈+ 2𝑟θ ̇) = -𝑞𝐵𝑟̇+ 𝐹 𝐸𝜃 (5.2) } Considering (5.
2) and rewriting it,

𝑚𝑟𝜃 ̈+ 2𝑚𝑟θ ̇+ 𝑞𝑟Ḃ(𝑡) = 𝐹 𝐸𝜃 And 𝑚𝑟𝜃 ̈+ 𝑟 𝑞𝐵 ̇(𝑡) 2 + 2𝑚𝑟θ ̇+ 𝑞𝑟Ḃ(𝑡) = 𝑟 𝑞𝐵 ̇(𝑡) 2 + 𝐹 𝐸𝜃 Consequently, 𝑟 2 (𝑚𝜃 ̈+ 𝑞𝐵 ̇(𝑡) 2 ) + 2𝑟𝑟̇(𝑚𝜃 ̇+ 𝑞𝐵(𝑡) 2 ) = 𝑟 2 𝑞𝐵 ̇(𝑡) 2 + 𝑟𝐹 𝐸𝜃
By integration with respect to the time, we obtain

𝑟 2 (𝑚𝜃 ̇+ 𝑞𝐵(𝑡) 2 ) = ∫(𝑟 2 𝑞 𝐵 ̇(𝑡) 2 + 𝑟𝐹 𝐸𝜃 ) 𝑑𝑡 (6)
We have now to solve this equation. As usual in Physics we begin to have a look on the limiting cases:

-If magnetic field is null, we obtain 𝜃 ̇𝑟2 = 𝐶 Which indicates a law of central force if 𝐶 ≠ 0. We choose thus 𝐶 = 0.

-If magnetic field is constant we obtain

𝑟 2 (𝑚𝜃 ̇+ 𝑞𝐵 2 ) = 𝐶
Consequently, it exists a system of coordinate where we can write the law

𝜃 ̇= - 𝑞𝐵 2𝑚 (7)
i.e. the frequency of Larmor [2].

-Considering now the case where the magnetic field is time depending. We add the relation (1) which becomes here

(∇ × 𝐸 ⃗ ). 𝑧 = 1 𝑟 [ 𝑑 𝑑𝑟 (𝑟𝐸 𝜃 ) - 𝑑 𝑑𝜃 𝐸 𝑟 ] = - 𝑑 𝑑𝑡 𝐵(𝑡) Using (4.3) 𝑑 𝑑𝜃 𝐹 𝐸𝑟 = 0 (8) Wich leads to 𝑑 𝑑𝑟 (𝑟𝐸 𝜃 ) = -𝑟𝐵 ̇(𝑡) Thus 𝐸 𝜃 = - 𝑟𝐵 ̇(𝑡) 2 (9)
And we obtain the system

{ 𝑟 2 (𝑚𝜃 ̇+ 𝑞𝐵(𝑡) 2 ) = ∫(𝑟 2 𝑞 𝐵 ̇(𝑡) 2 + 𝑟𝐹 𝜃 ) 𝑑𝑡 𝐹 𝐸𝜃 = -𝑞 𝑟𝐵 ̇(𝑡) 2 } Thus 𝜃 ̇= - 𝑞𝐵(𝑡) 2𝑚 ( 10
)
It is this relation which allows to build the model, without adding particular assumption.

Verification

Using the laws of classical mechanics, we verify if the solutions we obtained are coherent. A method is to use the power theorem. Indeed, we know that magnetic force doesn't 𝐸 𝜃 = 𝑚𝑟 𝜃 q (12.2) }

Potential vector

As presented we use the system given with (2) , (3) and (4.2). We rewrite it

{ 𝐸 ⃗ = - 𝛿 𝛿𝑡 𝐴 ∇ × 𝐴 =
{ ∇ × 𝐴 = 𝐵(𝑡) ⃗⃗⃗⃗⃗⃗⃗⃗ ∇ × 𝐸 ⃗ = - 𝑑 𝑑𝑡 𝐵 ⃗ 𝐸 ⃗ = - 𝛿 𝛿𝑡 𝐴 }
And with the 2 nde law of Dynamics. The force of Lorentz can thus here be rewritten

𝐹 = 𝑞𝐸 ⃗ + 𝑞𝑉 ⃗ × 𝐵 ⃗ = -𝑞𝑟 𝐵 2 𝑒 𝜃 ⃗⃗⃗⃗ + 𝑞𝑉 ⃗ × 𝐵 ⃗ (17)

Equations of motion

They are given by equations of Dynamics (5.1), which becomes

𝑚(𝑟̈-𝑟𝜃 ̇2) = 𝑞𝐵𝑟𝜃 ̇
And, with (10)

𝑟̈= - 𝑞 2 𝐵 2 4𝑚 2 𝑟
Thus a differential equation depending on the time

𝑑 2 𝑑𝑡 2 𝑟(𝑡) = - 𝑞 2 𝐵 2 (𝑡) 4𝑚 2 𝑟(𝑡) ( 18
)
Considering firstly the limiting case where the magnetic field is constant, we obtain an equation under the form

𝑑 2 𝑑𝑡 2 𝑟(𝑡) = -𝐶𝑟(𝑡)
Where 𝐶 is constant and whose solution is

𝑟(𝑡) = 𝐶 1 sin(√𝐶𝑡) + 𝐶 2 𝑐𝑜𝑠(√𝐶𝑡)
But 𝑡 and 𝜃 are here linked with a relation of proportionality, such

𝑟(𝑡) = 𝐶 1 sin(𝐶′𝜃) + 𝐶 2 𝑐𝑜𝑠(𝐶′𝜃)

Trajectory is thus a circle which passes through the origin of the system of coordinate, for example And speed of the point particle is constant. We reobtain thus a classical result of electrodynamics, but in other manner, and as a limiting case. We can now use the model for time varying magnetic fields.

-Simulations

As we introduced it we consider the magnetic field

𝐵 = 𝐵 0 + 𝐵 1 𝑒 𝑗𝑤𝑡𝑡 Thus 𝜃 = - 𝑞 2𝑚 [𝐵 0 𝑡 + 𝐵 1 𝑗𝑤 𝑒 𝑗𝑤𝑡 ] ( 19
)
Equation ( 18) becomes

𝑑 2 𝑑𝑡 2 𝑟(𝑡) = - 𝑞 2 [𝐵 0 + 𝐵 1 𝑒 𝑗𝑤𝑡 ] 2 4𝑚 2 𝑟(𝑡)
A solution is obtained using the functions of Whittaker

𝑟(𝑡) = 𝐶 1 𝑊ℎ𝑖𝑡𝑡𝑎𝑘𝑒𝑟𝑀 [- 𝛼 𝑤 , 𝛼 𝑤 , 𝛽] + 𝐶 2 𝑊ℎ𝑖𝑡𝑡𝑎𝑘𝑒𝑟𝑊 [- 𝛼 𝑤 , 𝛼 𝑤 , 𝛽] √𝑒 𝑗𝑤𝑡 with { 𝛼 = 1 2 𝑞𝐵 0 𝑚 𝛽 = 𝑞𝐵 1 𝑒 𝑗𝑤𝑡 𝑤𝑚 }
Where 𝛼 is the frequency of Larmor corresponding on 𝐵 0 . We study two cases:

6.1 𝐵 0 = 0
We obtain here

{ 𝜃 = - 𝑞𝑒 𝑗𝑤𝑡 2𝑗𝑚𝑤 𝛼 = 0 𝛽 = -2𝑗𝐵 1 𝜃}
To do a simulation we choose constants as follows : 𝑚 = 𝑞 = 𝑤 = 𝐵 1 = 𝐶 1 = 1 and 𝐶 2 = 0

Trajectory is this time a spiral, which passes again through the origin of the system of coordinate.

Considering now the general case. Solving (19) we obtain the expression of the time

𝑡 = - 𝜃 𝛼 + 𝑗 𝑤 𝐿𝑎𝑚𝑏𝑒𝑟𝑡𝑊 [ 𝐵 1 𝐵 0 𝑒 -𝑗 𝜃𝑤 𝛼 ⁄ ]
We obtain a family of curves, that we can again draw in a polar system of coordinate. To present briefly this family, we can note that the general form is linked to the ratio

𝑛 = 𝑤 𝛼
When this ratio is integer, the corresponding trajectories are closed, which isn't necessary the case for arbitrary values. Consider for example the following real parts of the solutions, drawn with

𝑚 = 𝑞 = 𝐵 1 = 𝐵 0 = 𝐶 1 = 𝐶 2 = 1 𝑅𝑒(𝑟(𝜃)) for 𝑛 = 1 𝑛 = 2 𝑛 = 3 𝑛 = 4 𝑛 = 5 𝑛 = 6
So we obtain a remarkable family of curves, which seems not have been studied before, and which present, as visible, interesting geometrical properties. Note their exhaustive and mathematical study won't be done in this paper, but we can present certain interesting imaginaries parts of the curves, such, for example, 𝐼𝑚(𝑟(𝜃)) for 𝑛 = 4 𝑛 = 6

-The mono dimensional harmonic oscillator case

Mono dimensional harmonic oscillator describes a straight. For this reason constant of area 𝐶 evoked in part 3 of the paper is again null. Consequently, the corresponding equation of motion is simply given with

𝑑 2 𝑑𝑡 2 𝑟(𝑡) = -𝑘𝑟(𝑡) - 𝑞 2 𝐵 2 (𝑡) 4𝑚 2 𝑟(𝑡)
Where 𝑘 is the coefficient of the attractive restoring force. As previously, we solve it in the case where magnetic field is constant. Solution is

-Discussion

Two points seems in favor of the validity of the model: Firstly, the good agreement between the expected behavior of the charged particle and the external magnetic field. Secondly, the fact the limiting cases can be reobtained using the equations of motion. However, if proximate works have already been published, [START_REF] De Jesus | Classical and quantum mechanics of a charged particle in oscillating electric and magnetic fields[END_REF], it seems difficult to compare results. Indeed, we notice important differences between our model and these works. One reason is they generally chose a potential vector ( [START_REF] De Jesus | Classical and quantum mechanics of a charged particle in oscillating electric and magnetic fields[END_REF], [START_REF] Kosmas | Charged Particle in an Electromagnetic Field using variational integrators[END_REF]), given with

𝐴 = - 1 2 𝑟 × 𝐵 ⃗
Because it satisfies relation (3). But, having now a look on our own expression, this time obtained from a demonstration, and given with ( 16)

𝐴 = - 1 2 [𝑟𝑒 𝑟 ⃗⃗⃗ + 𝑑𝑟 𝑑𝜃 𝑒 𝜃 ⃗⃗⃗⃗ ] × 𝐵 ⃗
We note the addition of a term between the expressions. Consequence is naturally a modification of equations of motion.

Future work

The question of the extension of the study of trajectories and their physical properties (speeds, energies, etc..) will naturally be an important way of possible later works. It seems us especially interesting to explore the phenomenon of resonance we have noted between the frequencies of Larmor and of variable magnetic field, and to predict possible measurable effects, able to be used in applied Physics. A purely mathematical approach can be to study the properties of the curves we obtained. At end, another works can extend the results to a quantum approach, using expression (16) to obtain the Lagrangian of the system.

Summary and conclusion

We have presented a model to describe the trajectories of a charged particle in an oscillating magnetic field, from the point of view of classical physics and electromagnetism. The equations of motions lead to original trajectories, which are closed when the ratio between frequencies of Larmor and of magnetic field is a rational number. Furthers works should compare these results with experimental data, and extend them to others theoretical approaches, in particular to quantum mechanics.

)

  Which is corresponding on the family of Roses curves, (i.e. a part of the Centered Trochoïde curves [3] as we obtained it elsewhere in any other way [4]. See for example a possible solution (𝑘 = 3.75).Consider now the general case. Solving equation of motion, part of this solution leads again to an original family of curves, depending on the choice of coefficients. See for example, for 𝐵 1 = 0.5 and 𝑘 = 1.25, the simulation


Marie-Pierre Pacaux-Lemoine

Quentin Berdal

Simon Enjalbert

Damien Trentesaux
email: damien.trentesaux@univ-valenciennes.fr

Towards Human-based Industrial Cyber-Physical Systems

Keywords: human-machine cooperation, adaptive level of automation, collaborative robotics, intelligent manufacturing systems, multi-agent systems, cyber-physical systems

The constant advances in sciences and technologies encourage industrialist and researchers in manufacturing, to address new challenges relevant to industrial Cyber-Physical Systems (CPS). Human aspects, among others, are of importance and researchers try to take them into account, but they remain to be efficiently dealt with during the design of industrial CPS. The goal of this paper is to highlight how it is possible to integrate "humanin-the-loop" inside the process control of industrial CPS. For that purpose, relevant studies already done in industrial engineering and Human-Machine Systems are presented, completed by an overview of the main cognitive dimensions industrial designers have to integrate in assistance systems definition, in order to benefit from human competencies and capacities while respecting human limits. The main idea is to balance Human and technology involvement, taking advantage of industrial CPS advances and Human capabilities identified and implemented through Human-Machine

Cooperation principles. The project HUMANISM, which is presented, aims to specify and experiment such principles.

I. INTRODUCTION

The constant advances in sciences and technologies encourage industrialist and researchers in manufacturing, to address new challenges relevant to industrial Cyber-Physical Systems (CPS). Human aspects, among others, are of importance and researchers try to take them into account. For example, Human Cyber-Physical systems have been defined as well as the notion of Operator 4.0 [START_REF] Romero | The Operator 4.0: Human Cyber-Physical Systems & Adaptive Automation Towards Human-Automation Symbiosis Work Systems[END_REF]. Control 5.0 has been also introduced, bringing together Control 1.0 to Control 4.0 in the artificial world, supporting the virtual-real duality as well as the Cyber-Social-Physical Spaces notion [START_REF] Wang | Control 5.0: From Newton to Merton in Popper's Cyber-Social-Physical Spaces[END_REF]. Most of the time, the design of an industrial CPS is mainly focused on how providing information to Human operators [START_REF] Lee | A cyber-physical systems architecture for Industry 4.0-based manufacturing systems[END_REF]. More, human aspects are often taken into consideration only at the end of the design process, that is once the control system has been fully designed. In fact, researchers in manufacturing mainly adopt a technocentered approach, favoring in priority the definition and allocation of tasks to automated intelligent systems while considering at the same time that humans, being operators or supervisors, will be there only to handle any unexpected situations efficiently.

Human aspects remain to be efficiently dealt with during the design of industrial CPS: in normal conditions, human operators or supervisors are assumed to decide and inform perfectly within due dates, with no stress or mental overload. They are called when a problem appears. They have to make quick and efficient decisions while they are facing complex systems composed by autonomous entities connected each other, with possible emergent behavior. Human operators thus face difficult tasks without being aware of process design or current state [START_REF] Pacaux-Lemoine | Designing Intelligent Manufacturing Systems through Human-Machine Cooperation Principles: A Human-Centered Approach[END_REF]. Consequently, in a single word, these humans are assumed by designers to behave like "magic humans" [START_REF] Trentesaux | A human-centred design to break the myth of the ''Magic Human" in intelligent manufacturing systems[END_REF], which is neither realistic nor sustainable in the context of the future industrial CPS interoperating with humans.

The goal of this paper is to highlight how it is possible to integrate "human-in-the-loop" inside the process control of industrial CPS. For that purpose, relevant studies already done in industrial engineering are presented in the next part, completed by an overview of the main cognitive dimensions industrial designers have to integrate in assistance systems definition, in order to benefit from human competencies and capacities while respecting human limits. The third section presents the Human-Machine System (HMS) approach, which deals with the support of the cooperation between human operators and assistance systems. An example of the use of HMS in the robotics field is detailed not only to illustrate this approach but also to identify its possible benefits in the context of human-based industrial CPS. To illustrate how human-based industrial CPS can be designed, the fourth section presents the project HUMANISM, which aims at implementing HMS approach to three types of CPS taking part in a Flexible Manufacturing System (FMS).

II. ISSUES RELEVANT TO HUMAN ASPECT IN INDUSTRIAL

CPS A growing number of researchers, especially in ergonomics and human-engineering, have already addressed the domain of industrial engineering to ensure more human-centered designs of manufacturing control systems. As a prerequisite, they have worked on Levels of Automation (LoA) and Human-Automation Symbiosis [START_REF] Lindström | Aligning manufacturing strategy and levels of automation: A case study[END_REF] [START_REF] Romero | Towards a Human-Centred Reference Architecture for Next Generation Balanced Automation Systems: Human-Automation Symbiosis[END_REF]. They propose methods and associated assistance systems in order to identify and design appropriate LoA, regarding human competences and needs [START_REF] Oborski | Man-machine interactions in advanced manufacturing systems[END_REF][8], as well as industrial system goals. Targeted industrial studies involving human aspects aim at providing cost-efficient solutions based on flexibility and functionality of proactive assembly systems [START_REF] Dencker | Proactive assembly systems-realising the potential of human collaboration with automation[END_REF], agile and adaptive manufacturing control architecture [START_REF] Leitao | ADACOR: A holonic architecture for agile and adaptive manufacturing control[END_REF], and they aim at integrated physical/mechanical and cognitive/information-related variables [START_REF] Frohm | Levels of Automation in Manufacturing[END_REF]. Other interesting studies have been conducted focusing on human-machine interface like augmented reality technics [START_REF] Malý | Augmented reality experiments with industrial robot in industry 4.0 environment[END_REF], or monitoring for maintenance [START_REF] Galar | Context awareness for maintenance decision making: A diagnosis and prognosis approach[END_REF], or are focusing on emerging technologies such as cobots [START_REF] Cherubini | Collaborative manufacturing with physical human-robot interaction[END_REF]. Applications to manufacturing cells are also found [START_REF] Pires | Semi-autonomous manufacturing systems: The role of the human-machine interface software and of the manufacturing tracking software[END_REF] [START_REF] Solvang | Shop-Floor Architecture for Effective Human-Machine and Inter-Machine Interaction[END_REF].

Industrial CPSs are dynamic systems regarding the different constraints they have to respect and especially temporal pressure, reaction time lags and numerous variables that could lead to unstable behaviors. The use of an HMS tools and approaches may then counter-balance the techno-centered design way followed by industrial CPS engineers. Among these tools and approaches, a first interesting one is the use of the five cognitive dimensions influencing human operator in the control and the supervision of dynamic situations [START_REF] Hoc | Cognitive Evaluation of Human-Human and Human-Machine Cooperation Modes in Air Traffic Control[END_REF]. These dimensions are the supervisory span, the control directness, the process information accessibility, process speed and continuity.

The first dimension relates to the supervisory span. It concerns the ability for human to access to process variables due to temporal, causal or spatial restriction or excess and their combination. A restricted supervisory span may result in the difficulty to anticipate process states regarding the speed of action feedbacks. An example of a dynamic situation presenting this characteristic is the control of heavy and long ship. A ship needs time to reach the direction or the speed defined by the pilot, and external unpredictable events like waves, winds and sea currents disturb its trajectory. Assistance systems provide useful support for anticipating disturbances using simulation tools. On the opposite, large supervisory span may lead to difficulty to identify process, to build a model and to make decision at the right time before decision begins to be obsolete. Typical example is fighter aircraft piloting. During certain phases of the flight, a pilot is unable to make a decision due to the quickness of the aircraft and the quick sequences of actions. All actions must be planned before the flight and updated during the flight. Recent fighter aircrafts are now equipped with assistance systems aiming at supporting task plan update on the tactical situation (SITAC display). This dimension deals with system complexity regarding the difficulty to identify causal relationship between variables. Industrial CPS and especially Intelligent Manufacturing Systems integrating smart products would present such a large supervisory span. Indeed, smart products have a quick and reactive behavior that impede human operator to intervene correctly for updating parameters. Assistance system based on simulation allows the test of commands regarding the current state of the smart products, and/or can provide a selection of possible commands to the human operators [START_REF] Pacaux-Lemoine | Designing Intelligent Manufacturing Systems through Human-Machine Cooperation Principles: A Human-Centered Approach[END_REF]. Another type of example is the lack of process awareness, when different staffs, like maintenance, designer and production staffs, modify the system regarding different and conflicting objectives. In that case, process may have an unexpected behavior because the model used by the human operator is not the good one anymore. Supports about the situation awareness using collaborative tools may mitigate such risks [START_REF] Millot | A common Work Space for a mutual enrichment of Human-machine Cooperation and Team-Situation Awareness[END_REF].

The second dimension relates to the control directness. It deals with the length of causal chains and the impact of intermediate variables, not fully controlled or predicted. Control directness affects system controllability, mainly due to a large number of variables. A large control directness, with too long casual chains, may prevent human operators from anticipating what could be the impact of their actions. A process with restricted control is more robust because disturbances have less chance to appear. Control directness relates then to complicated and dynamic systems. Industrial CPSs can be regarded as these kind of systems since they are concerned with numerous intermediate variables, with several communications between entities with sometimes hidden connections or connections unknown by human operators. Planes, trains, plants present large control directness, and more and more assistance systems have to be designed in order to restrain this dimension to allow human operator to grasp global process through its sub-parts. HMS provides tools to make shorter causal chains with assistance systems like information screening, prediction and simulation supports.

The third dimension concerns the process information accessibility. This dimension also concerns accessibility of process variables but it deals with the absence of the variable that must be calculated or assessed in the base of "surface" variables. System observability is relevant to this dimension. A variable, which is not accessible, implies that the human operator has to make inferences and hypotheses. It could also be a problem of delay. Indeed, a variable can be accessible but too late to diagnose process and make decision. In this case, predictions based on model have to be done. One of the advantages of industrial CPS concerns the potential access and use of large amount of raw data provided by sensors connected to physical components of the process (big-data). This implies then that diagnosis assistance systems are provided to support maintenance and supervision decisions.

The fourth dimension relates to process speed. It concerns typically sampling frequency and the time human operator has to monitor processes to be sure not to miss an important information. This dimension also deals with the human operator's ability to plan during control task if the process is slow, e.g. blast furnace, or fast like in highly automated systems. In the first case, simulators provide useful tools to allow human operator to build a representation of the process according to a specific time. In the second case, assistance systems like autopilot replaces human operator's skills and reduces his/her workload due to time pressure induced by process speed.

The last dimension concerns process continuity that is the evolution of the process and the ability for human operator to determine significant variations. Assistance systems aim at highlighting important information stemming from data analysis, especially regarding sampling.

Human operators feel these five dimensions differently according to their own expertise, experiences and capacity, and this last ability can evolve very quickly according to the difficulty of the tasks.

Human-Machine Systems approach uses these five dimensions in order to design assistance systems adapted to process characteristics and human operators' needs. Such an approach would be very useful in the design of Human based industrial CPS. A detailed presentation of this approach is now provided in the next section.

III. THE HUMAN-MACHINE SYSTEMS APPROACH

Industrial systems designers have now to work with computer sciences and psychological designers, and they are at the first steps of identifying their common objectives, their sharing languages and models. The construction of cooperation between researchers who belong from different fields is a long way which have already been followed by the Human-Machine Systems (HMS) domain. HMS approach aims at considering Human operator and assistance system in the same way, each one having own competences and capacities which have to be update if necessary to enable fruitful cooperation. HMS first principles have been proposed in 1975 by Thomas Sheridan who described models of human performance concerning information, control and decision [START_REF] Sheridan | Man-Machine Systems: Information, Control, and Decision Models of Human Performance[END_REF]. Such principles have the objectives to adopt a Human-Centered design approach, considering the current and future activity of Human operator at the first steps of the HMS design. Human-Machine Cooperation focuses on cooperative goals and aims at taking into account objectives, competences and characteristics of each agent, human and artificial ones, to find the best organization regarding the situation to control [START_REF] Hoc | Cognitive Evaluation of Human-Human and Human-Machine Cooperation Modes in Air Traffic Control[END_REF].

A. Human-Machine Cooperation principles

The principles of Human-Machine Cooperation have been defined more than thirty years ago [START_REF] Millot | An approach of dynamic allocation of supervision tasks between man and computer in control rooms of automatized production systems[END_REF] and been used in different domains like nuclear plant, air traffic control, car driving, military robotics and aviation [START_REF] Millot | A common Work Space for a mutual enrichment of Human-machine Cooperation and Team-Situation Awareness[END_REF] [START_REF] Flemisch | Shared and cooperative control of ground and air vehicles: Introduction and General overview[END_REF]. The general purpose is to propose a generic method to assist human designers in identifying how the human operator and machine may interact with suitable adaptive levels of automation while guaranteeing their performance, safety and security [START_REF] Habib | A method for designing levels of automation based on a human-machine cooperation model[END_REF]. Human-Machine Cooperation studies led first to define a model of a cooperative "agent", a human being or an artificial entity, with two main dimensions: -the agent's ability to control the process, also called the Know-How (KH), and -the agent's ability to cooperate with other agents involved in the process control, also called the Know-How-to-Cooperate (KHC).

The Know-How of an agent only concerns the control of the process, the achievement of individual tasks without taking into account potential interaction with other agents. The KH is split up into two parts, one is called the internal KH, the other the external KH. The internal KH relates agents' competences and capacity to control the process. The competence of a Human agent is mainly composed by knowledge, rules and skills to control the process [START_REF] Rasmussen | Skills, Rules, and Knowledge; signals, Signs, and Symbols, and Other Distinctions in Human Performance Models[END_REF]. It is linked to expertise, experience and practices of agents. Similarly, the competence of an artificial agent is mainly related to its ability to acquire knowledge (like with machine learning), to follow rules (like expert systems) and to apply predefined commands (like with PLC, Petri net). The capacity of a Human agent is mainly related to workload, fatigue and attention, while an artificial agent is constrained by energy, memory and processor capacity. The external KH deals with the ability to get information from the process and the ability to act on the process. Therefore, the internal KH seems to be close to the cyber part of an industrial CPS and the external part seems to be close to the physical part of the industrial CPS. For implementation purpose, the KH has been summed up in four functions: information gathering, information analysis, decision making and action [START_REF] Parasuraman | A Model for Types and Levels of Human Interaction with Automation[END_REF], and these functions only concern the process and not how the activity of other agents involved in the process control can be considered. The KHC is dedicated to this special task.

The Know-How-to-Cooperate is also split up into two parts. The external KHC is the ability of an agent to have information about other agents and to provide information to other agents. The support of the external KHC is called the Common Work Space. It supports the situation awareness dedicated to process state and environment, but it is enriched by the team situation awareness dealing with past, current and future activity of all agents [START_REF] Millot | A common Work Space for a mutual enrichment of Human-machine Cooperation and Team-Situation Awareness[END_REF]. The internal KHC allows agents to build up a model of others in order to make easier the cooperation with them. It is built up and updated by learning, training and exchanging with agents. Agents gather and analyze information about others in order to infer their KH and KHC. This notion seems to be inexistent in industrial CPS where other artificial agents are most of time considered as part of the process or they have myopic behavior combined to a reactive behavior that impedes Human operators to have deep exchanges with them.

KH and KHC can be understood as two parallel functions described each one by different sub-functions. Fig. 1 contains an illustration of this articulation. In this figure, the cooperative agent model is used to highlight interaction between Human and machine [START_REF] Pacaux-Lemoine | Towards vertical and horizontal extension of shared control concept[END_REF]. KH functions of agents are in interaction by the means of a Common Work Space represented by the blue area in the middle of the figure. KHC functions of agents use information provided by the Common Work Space in order to build up a model of the other agent and to evaluate this agent's involvement in the process control. The results of the evaluation compared to the own involvement and model of oneself (interference detection and management) have to lead to the adjustment of the position of the four sliders that describe the functions allocation (scales represented in the Common Work Space). A scale is associated to each KH functions. The position of the slider on the scale defines the degree or the percentage of sharing between Human and machine activity. Some functions can be completely allocated to machine or to human regarding their competency and capacity, but also regarding how they can take into account the activity of the other. The function allocation can be predefined and updated according to the information on the current situation. Model of cooperative activities (adpated from [START_REF] Pacaux-Lemoine | Towards vertical and horizontal extension of shared control concept[END_REF]).

The model is also able to encompass existing definitions of levels of automation (LoA) and to propose new ones [START_REF] Habib | Adaptation of the level of automation according to the type of cooperative partner[END_REF]. It is generic enough to be used on different types of levels of activity, also called layers [START_REF] Pacaux-Lemoine | Layers of Shared and Cooperative Control, assistance and automation[END_REF]. Layers are most of the time used at operational, tactical and strategic levels. Industrial CPS are usually more concerned by layers close to command/control of process, than higher decisional layers close to the strategic level, where decision making is more and more distant from the physical part of the process.

B. An illustrative example

To illustrate these principles and to visualize their potential benefits to design industrial CPS that integrates efficiently the human operator or supervisor, an application dealing with crisis management and cooperation between a supervisor and robots is presented. This application, designed and realized in authors' lab, presents some similarities with the context of industrial CPS and research results can thus be extrapolated. For example, in industrial CPS, one can easily imagine that several human operators augmented with information capabilities (smart glasses, smart garment) share spaces and activities with mobile autonomous robots or intelligent products. In this context, crisis could concern hazardous situation, critical machine breakdown management, threat or unexpected process states that appear.

In this study, three types of agents take place in the scenario: a human supervisor (equipped with EEG and a brain computer interface BCI), operators (equipped with smart garment) and robots [START_REF] Pacaux-Lemoine | Human-Robots Cooperation through Intelligent Garment[END_REF]. At the operational level, robots and operators have to be cooperative with each other and with the human supervisor in order to achieve a given common goal (solving the crisis or limiting its impact). The human supervisor is at the tactical level; he/she can observe the situation, make decision about the process and send instructions to the robots using a visual support presenting large amounts of information [START_REF] Habib | Adaptation of the level of automation according to the type of cooperative partner[END_REF]. The model of cooperative activities helps agents to better organize their actions. Fig. 2 proposes an example of the use of the model of cooperation for identifying and selecting shared control functions and common work space in order to support cooperation between the supervisor and one robot in the field [START_REF] Pacaux-Lemoine | Human-Robot Cooperation through Brain-Computer Interaction and Emulated Haptic Supports[END_REF]. Four LoA have been identified according to KH and KHC of each agent (cf. Fig. 2). The human KH is violet, the robot KH is blue, the human KHC is light pink and the robot KHC is light blue. If the agent has the ability to perform one of the KH or KHC functions presented in the model (cf. Fig. 1), "1" is written in the rectangle. A new concept, the "emulated shared control" has been proposed to improve the robot KHC. This new concept has been defined with the same idea as the haptic control used in robotics and aeronautics (force feedback in the joystick) and in car driving (force feedback in the steering wheel) [START_REF] Abbink | Haptic shared control: smoothly shifting control authority?[END_REF]. In the case of a command sent via the BCI, there is no real haptic feedback because no muscles are involved in the control to oppose or to follow the direction provided by the system. But the idea is to emulate this haptic behavior with a visual display. Depending on where obstacles are detected, the robot/BCI system makes it easier or more difficult-in terms of (mental) effort-for the human to deliver a command (cf. Fig. 2 where blue and red arrows give the left or right direction to the robot). The common work space focuses on this visual feedback and contains a live video stream of the environment from the robot's perspective. Red rectangles highlight the differences between each LoA regarding the layer, operational or tactical. A command is associated to each "1" function. Model predictive control approach (MPC) has been used to combine the different functions and their priority [START_REF] Habib | Adaptation of the level of automation according to the type of cooperative partner[END_REF]. Each LoA has been evaluated within the framework of an experiment conducted up to now with only one participant and one robot. Comparison between results recorded for each LoA highlights the interest of the LoA defined by "With Emulated Haptic Feedback" and "Without Obstacle Avoidance". Precise analysis of agents KH and KHC during each step of the experiment underlines lack of situation awareness regarding obstacle detection and avoidance. This result leads to improve agents KHC to share their situation awareness. Shared control and common work space (adapted from [START_REF] Pacaux-Lemoine | Human-Robot Cooperation through Brain-Computer Interaction and Emulated Haptic Supports[END_REF])

The robot used in this experiment owns very simple KH functions, and difficulty holds here in the interaction process and the design of its KHC. Obviously, the control of a swarm of artificial agents increase considerably operator's and supervisor's tasks. From the knowledge drawn from these experimentations, we concluded that the complexity to elaborate efficient and effective human-based industrial CPS will mainly come from the number and the variety of artificial agents involved. Indeed, in the context of industrial CPS, smart products, smart machines, smart inventories, smart tools, etc. would compose such artificial agents. From our perspective, the complexity of the HMS will be brought by the complexity of human tasks (e.g. supervision and control of a swarm of smart entities), by the inherent complexity of each of these smart entities (e.g. issuing their specific supervision) and complexity in the interaction between human and all the set of these artificial agents (e.g. occurrence of emerging behaviors and self-organization).

To study some aspects relevant to these kinds of complexity is the topic of the French ANR granted project HUMANISM, which is presented in the following section.

IV. THE HUMANISM PROJECT

The objective of the HUMANISM project, newly launched in October 2017, is to study, improve and experiment innovative and cooperative systems involving production and supervisory human operators, diagnostic systems and control systems of intelligent Flexible Manufacturing Systems (FMS) within the context of industrial CPS. Two complementary contributions are proposed and concerned operational and tactical levels. The first one deals with the design of a method to control the intelligent FMS assuming that it is composed of autonomous artificial agents (which are known to provide adaptation capacity and robustness [START_REF] Barbosa | Dynamic selforganization in holonic multi-agent manufacturing systems: the ADACOR evolution[END_REF]) through a better integration of the Human operator. Complementary and new approaches have to be proposed and developed to enable the cooperative diagnostic of technical incidents and the control of emergent behaviors of such systems. The second contribution focuses on the identification of the pertinent information set that allows the design of monitoring and control systems enabling to maintain human operator's situation awareness in order to help him/her to take (near)optimal decisions, and to define and select LoA that takes into account competences and capacities of each agent involved in the global system, being human or artificial.

The cooperation models and methodology presented in the previous section will be used, as well as an approach called the Cognitive Work Analysis (CWA). CWA is a method to support the analysis of decision-making process and of the evaluation of situation awareness, workload and modes of cognitive control in dynamic situations [START_REF] Chauvin | Using Cognitive Work Analysis to design smart grid interfaces[END_REF]. To ensure the genericity of our theoretical developments, HUMANISM will be applied on three different kinds of autonomous systems, called Artificial Self-Organized systems (ASO), with a focus on their ability to contribute to a controlled and efficient intelligent FMS: a cobot system, a swarm of intelligent products, and a swarm of mobile robots (cf. Fig. 3). They present different adaptable or adaptive LoA. They differ each other in the kind of interaction they have with the Human (close vs. distant), and in their behavior (predictive or reactive, with or without communication abilities and accessible data). Multi-level cooperation between Humans and innovative autonomous systems

From a theoretical point of view, HUMANISM will give the opportunity to observe and to analyze the propagations of constraints between the operational level and the tactical level of an intelligent socio-technical system as well as their impact on operators' activities.

Experiments will be conducted with several participants on the three ASOs combined to the support systems of the tactical level. Of course, experiments will be counter-balanced regarding learning and order effects. Subjective data will be stemmed from questionnaires, interviews and video records in order to analyze participants' activity. Objective data (duration of operation, number of interferences between human and ASOs, number of mistakes...) will reinforce the results in terms of situation awareness, acceptability, performance, safety.

V. CONCLUSION In this paper, authors fostered the idea to integrate the principles of Human-machines systems to improve the effectiveness of the cooperation between Human, being operators or supervisors, and industrial Cyber Physical Systems composed of autonomous artificial agents. Pending issues concerning integration of Human in industrial Cyber Physical Systems and associated cognitive dimensions have been discussed. The lack of tools facilitating this integration and especially the lack of principal axes [START_REF] Wittenberg | Human-CPS Interaction -Requirements and humanmachine interaction methods for the industry 4.0[END_REF], leads us to investigate principles of cooperation in Human Machine Systems. An illustrative and concrete example helped us to identify the different kinds of complexity that may make it difficult to design human-based industrial Cyber-Physical Systems that will be studied within the context of the HUMANISM project. Both approaches from Human-Machine Systems and industrial Cyber Physical Systems will be mixed and should be evaluated to assess performance of human-based industrial CPS.

Fig. 1 .

 1 Fig. 1.Model of cooperative activities (adpated from[START_REF] Pacaux-Lemoine | Towards vertical and horizontal extension of shared control concept[END_REF]).

Fig. 2 .

 2 Fig. 2.Shared control and common work space (adapted from[START_REF] Pacaux-Lemoine | Human-Robot Cooperation through Brain-Computer Interaction and Emulated Haptic Supports[END_REF])

Fig. 3 .

 3 Fig. 3.Multi-level cooperation between Humans and innovative autonomous systems

ACKNOWLEDGMENT

The research presented in this paper is carried out in the context of the HUMANISM N° ANR-17-CE10-0009 research program, funded by the ANR ''Agence Nationale de la Recherche", and by SUCRé project, partly funded by the Regional Council of the French Region "Nord-Pas de Calais-Picardie" and supported by GRAISYHM programs. The authors gratefully acknowledge these institutions.