1. Introduction

The electromagnetic modeling of the reflection by a given structure generally provides reflection coefficients which depend on several parameters: wavelength, polarization, observation angles. The color of the object also depends on the light source, and must take into account human color perception. In the first part of this paper, we present a convenient way to convert these reflection coefficients into colors that can easily be displayed on a computer screen, then printed.

In the second part of the paper, we apply this method to the study of the structural color of Morpho rhetenor butterflies. Their blue iridescent color is not due to pigments, but to an interference phenomenon involving a complicated grating structure. Previous studies in the same domain have approximated this structure using a stack of thin films [1,2]. Here, we take into account the microstructure of the wing, and model it as a multilayered grating of ribs. The problem is solved using a rigorous modal method which is both fast and accurate. It allows us to solve difficult problems like wings illuminated by skylight. Indeed, this problem requires the computation of the response for a huge collection of directions of incidence, for the two fundamental polarizations, and for a number of wavelengths inside the eye spectral bandwidth. All these data are then converted into colormaps representing the color observed for different observation directions. We also use our method to recover several experimental results obtained by Vukusic et al. [1] concerning the spatial and spectral distributions of the energy diffracted by the wing.

2. From reflection coefficient to color

In this section, our objective is to summarize the basic ideas that we have used for converting reflection coefficients into colormaps. Basically, we have picked out and put together several methods well known in color science and graphic imaging. We will therefore not give all the details, and this section should be considered as a collection of recipes. For a quick illustration of the main idea, we recommend the spectrum applet in Ref. [3].

Let us suppose that the body under study is illuminated by an illuminant characterized by its energy distribution D(λ). In this paper, we use the CIE (Commission Internationale de l–Éclairage) normalized illuminant D65, which closely matches that of the sky daylight [4–6]. For convenience, we have used a 5^th order polynomial approximation of this function (Fig.1a). Assuming that the body has a reflectivity R(λ), we can compute the CIE XYZ tristimulus values [7] which are one possible way to characterize the reflected color, according to:

where the weighting functions x̄, ȳ, z̄ (Fig.1b) are the CIE 1964 spectral tristimulus values [5,6,8] and k is a normalization factor defined in such a way that an object with a uniform reflectivity R(λ)=1 gives a luminance component Y equal to 1:

In the numerical treatment, all the integrals are replaced by discrete sums. Due to the human eye spectral bandwidth, the only values of the wavelength that we will retain lie in the range 0.38–0.78 µm. We use 11 equally spaced sampling values in this range, which are enough to get convergence in the color rendering.

 

Fig. 1. a: Relative energy distribution of the D65 illuminant (blue) and its 5th order polynomial approximation (red) versus the wavelength λ. b: Spectral tristimulus values.

Download Full Size | PDF

We then transform the XYZ components into RGB components, which are more convenient for the imaging process. Since RGB values are device dependent, we use the BT.709 Recommendation, which is representative of contemporary monitors in studio video, computing and computer graphics. In this case, the conversion is given by the linear transformation [7]:

These RGB values are normalized and should stay between 0 and 1. In fact, they can in some cases be slightly less than 0 or greater than 1. It means that the color can not be accurately represented on such monitor devices by additive reproduction of these RGB values. In this case, the values are rounded to 0 or 1.

The last step is the color representation. For this purpose, the PPM (Portable PixMap) file format [6] is very simple to use: it is a text file in which one simply enter the RGB values for each pixel. Note that each RGB value should be an integer between 0 and 255, obtained by scaling the values obtained from (3). These PPM files can be easily displayed on a computer screen using classical image editors, and easily printed.

At this stage, we have replaced a complicated set of data (i.e. the spectrum R(λ)) by a color point in a colormap. As illustrated in the following sections, this representation can be very useful for obtaining a global insight on the reflection properties of a body. A two-dimensional colormap based on this representation provides two free parameters (for example the two angles of diffraction).

3. Modeling the wing structure

Morpho rhetenor butterflies’ wings have a complicated structure composed of scales [9,10]. The cross-section of the ground scales (those producing color effects) exhibits a grating structure whose period is made of ridges with a tree-like geometry (Fig.2a). In the Morpho rhetenor butterfly scales, the lamellae (the “branches” of the tree-like structure) run near-parallel to the base of the scale [10]. So it is easier to model than some other species where the lamellae are tilted to the base of the scale [10,11], and this is the reason for which we focus on the Morpho rhetenor species. According to measurements reported in [1], the complex optical index of this tree-like structure will be approximated by the constant value n=1.56+i 0.06 in the whole optical range. Note that our model can take into account dispersion if the n(λ) data were available. The numerical analysis is performed with the help of a rigorous method dealing with stacks of lamellar grating layers. It is an x-invariant structure in which each grating layer is y-periodic and can contain several rectangular regions with different optical indexes (in the present study, there are only two regions with indexes n and 1, see Fig.2b). The structure is illuminated by a plane wave with arbitrary wave vector (Fig.3a). In each grating layer, the electromagnetic field is expanded on a modal basis [12–15], whereas in the superstrate and in the substrate (homogeneous media), the electromagnetic field is expanded on Rayleigh expansions (Fourier basis). On each interface between two gratings layers (red planes of Fig.2b), or between a grating layer and the substrate or the superstrate, the continuity of the electric and magnetic field components is expressed by projecting the fields expansions on the Fourier basis. Finally, one gets the characteristics of the field reflected and transmitted by the grating. In the present case, we are mainly interested in the reflected field. Each computation gives, for a given incident plane wave characterized by the parameters θ ⁱ , φ ⁱ , δⁱ , λ, the fraction of energy (the grating efficiency) reflected in the n^th order in the direction defined by ${\mathrm{\theta }}_n^r$ and ${\mathrm{\phi }}_n^r$ which is a function $e_{{\theta }_n^r,{\phi }_n^r}^r({\theta }^i,{\phi }^i,{\delta }^i,\lambda )$.

 

Fig. 2. a: Transmission electron microscope image showing the cross-section through a single Morpho rhetenor scale. Image reprinted from [1], with permission from P. Vukusic and the Royal Society. b: Modeled structure; the two red lines define a grating layer, the optical index of black regions is n, and that of white regions is 1.

Download Full Size | PDF

 

Fig. 3. a: k ⁱ is the incident wave vector, with angles θ ⁱ and φ ⁱ . s ⁱ and p ⁱ are respectively normal and parallel to the incidence plane Π ⁱ . The polarization of the incident electric field E ⁱ is defined by the angle δⁱ . b: Schematic representation of two conical subregions of the upper hemisphere.

Download Full Size | PDF

Let us focus on the following problem: we suppose that the butterfly wing is illuminated by the skylight, and we want to determine the color reflected by the wing in the direction (θ ^r ,φ ^r ). Of course, it is a particular problem and, using the same technique, we could also model the color reflected by the wing if it was for instance illuminated by a focussed beam. Since the wing is opaque and supposed to lie in the horizontal plane, we only need to consider the incident light coming from the upper halfspace. We divide this halfspace in conical subregions centered around a direction (θ,φ), with solid angle ΔΩ_θ,φ (Fig.3b). In order to model the uniform luminance of the skylight, we consider that two plane waves impinge on the wing with a wave vector k ⁱ colinear to the axis of each cone. These two plane waves have orthogonal polarizations (δⁱ =0 and 90°) since daylight is unpolarized. Depending on θ ⁱ , φ ⁱ and the grating period, they produce one or several reflected orders flowing inside one or several solid angles $\Delta {\Omega }_{{\theta }_n^r,{\phi }_n^r}$,. For each wavelength λ, and each conical subregion centered around a reflected direction (θ ^r ,φ ^r ), we add incoherently the contributions of all the incident waves giving a reflected order in this subregion, regardless of polarization, and obtain the power reflected in this subregion, which we denote by $R_{{\theta }^r,{\phi }^r}\left(\lambda \right)$. Using the recipe described in section 2 where R(λ) is now replaced by $R_{{\theta }^r,{\phi }^r}\left(\lambda \right)$, we obtain a color (including the brightness information) which can be denoted by a function RGB(θ ^r ,φ ^r ), and which can be easily plotted.

Note for instance that if the halfspace is divided in conical subregions with 10° aperture (i.e. 9 values of θ ⁱ and 36 values of φ ⁱ , with 2 values of δⁱ and 11 values of λ), one has to compute 9×36×2×11=7128 grating problems). Our computer code, running on a Pentium IV 1.4 GHz computer, deals this problem in about 18 minutes.

4. Numerical results

We first need to obtain convenient values for the various geometrical dimensions involved in the Morpho rhetenor scales. For this purpose, we have used the image of Fig.2a, and determined the scaling factor in order to match some measured values reported in Ref. [1]. The grating period is determined in order to have the same diffraction directions as in Fig.4a of [1], the scales being surrounded by air as well as by isopropyl alcohol (IPA), and illuminated in normal incidence with λ=488 nm. Since the scale structure differs from a real periodic grating, some diffusion is added and enlarges the diffraction peaks, which makes the determination of the grating period inaccurate. A more accurate value is obtained by matching the spectral behavior of the scale reflectivity (Fig.5a of [1]). Combining these two fits, we obtain a grating period of 746 nm, which will be used throughout this section. The other characteristic dimensions are shown on Fig.4.

 

Fig. 4. Characteristic dimensions of the modeled scale (in nanometers). The white regions are filled with air (or with a liquid solvent), and the gray regions are filled with a material with optical index n. In this figure, and using the tree image, we will say that each tree has four branches. The lower layer models the bottom membrane of the scale. Apart from the bottom air layer, all the dimensions are scaled.

Download Full Size | PDF

Using a structure with the same dimensions as those given in Fig.4, but with 16 branches for each tree, we compute the diffraction efficiencies (which are the fractions of the incident energy diffracted in the various grating orders) when the scale is illuminated in normal incidence by a plane wave with λ=488 nm. We take an average value of the efficiencies obtained for the two fundamental polarizations in order to simulate an unpolarized incident light. The computations are made when the scale is in air and when it is immersed in IPA (with optical index 1.38) and are plotted on Fig.5. When the scale is in air, the main part of the diffracted light is reflected in the ±1 orders, whereas when the scale is in IPA, the index contrast between the scale and the external medium is weaker, and almost all the light is transmitted in the 0 order. Note that the missing part of the energy is absorbed by the structure.

 

Fig. 5. Diffracted efficiencies by a scale in air and in IPA. The reflected angles are between 0 and 180° and the transmitted angles are between 180 and 360°. The normals to the scale correspond to 90° and 270°.

Download Full Size | PDF

Figure 6 shows the total reflection and transmission of the same structure (16 branches, lying in air) illuminated by a plane wave in normal incidence. These total reflection and transmission coefficients are obtained by adding all the reflected or transmitted efficiencies, for the two fundamental polarizations. We call E‖ (resp. H‖) the case where the electric (resp. magnetic) field is parallel to the scale ridges (x-axis in Fig.2b). Note that these cases are respectively denoted by TM and TE in Ref.[1]. The reflection curves are very similar to the experimental data reported in the Fig.5 of Ref.[1], and show a pronounced peak of reflectivity for the smaller wavelengths, giving to the wing its blue color. Note that the transmission curves predict that a sizeable amount of energy is transmitted for larger wavelengths, which does not agree with the experience. This phenomenon is quite general in our modeling: we generally find a good agreement where the amount of energy is significant, and a poor agreement where the amount of energy is low. This discrepancy can be due to several reasons, and we give below those which appear to be the most relevant.

First, we can observe that when the amount of energy is significant, the phenomenon is well described by a periodic structure and its interference process, whereas when the amount of energy is low, the diffusion by the random characteristics of the structure (which is not taken into account here) seems to give a significant contribution.

Secondly, the wing structure in the region between the ridges and the bottom membrane is actually rather complex and disordered, with probable inclusions of pigments. It will be shown below (Fig. 13) that for shorter wavelengths the field does not go deeply inside the structure, contrary to longer wavelengths. It means that the longer wavelengths are much more sensitive to the bottom region of the scale. If the absorption phenomenon in the actual wing is mainly due to localized effects inside the pigments (or other absorbing defaults) located in this bottom region, it explains the discrepancy between our computations and the experimental data for longer wavelengths. Note that in our model, the absorption is only due to the imaginary part of the ridges’ permittivity.

Nevertheless, since the color is generated by the stronger reflected components, we can be confident in the relevance of our approach as regards the color rendering.

 

Fig. 6. Total reflection and transmission for both polarizations.

Download Full Size | PDF

In order to present the color appearance of the scale, we now implement the approach exposed in section 2. We suppose that the wing is illuminated by a wide source modeling the sky, with the spectral energy distribution of the D65 illuminant, and assume that the incident light is coming from the whole upper halfspace. We divide this halfspace in conical subregions centered around directions (θ,φ) separated with angular variations of 10° degrees, then we collect and add the energies of the diffracted orders reflected in each conical subregion centered around a reflected direction (θ ^r ,φ ^r ). Finally, we transform this information in a RGB(θ ^r ,φ ^r ) map.

Figure 7 shows this map for the structure described previously (16 branches, lying in air). On these maps, the horizontal axis gives the θ ^r value, from 0 to 80°, and the vertical axis gives the φ ^r value, from 0 to 180°. Note that the definition of θ ^r and φ ^r follows that of Fig.3a. It means that the left side of the maps gives the reflected color in the direction normal to the scale (θ ^r =0) and the right side gives the reflected color for near grazing angles.

In figure 7, the left map color scale intensity goes from black (nothing reflected) to almost white (in fact the D65 illuminant color) if all wavelengths were totally reflected. This map appears dark since the reflection is generally much lower than unity. Since the human eye can distinguish between much more different levels of illumination than a computer monitor or a printer can render, this dark map do not really agree with the visual perception. In the right map, the global brightness is enhanced in such a way that the greater R, G or B value among all the map reaches the maximum value of 255. This map retains the relative brightness information between square subregions and agrees better with the visual perception.

 

Fig. 7. Reflected colors by a scale with 16 branches in air. The horizontal axis gives the θ ^r value, from 0 to 80°, and the vertical axis gives the φ ^r value, from 0 to 180°. The angular shift between two neighboring squares is equal to 10°.

Download Full Size | PDF

When the scale is immersed in IPA, the index contrast between the scale and the external medium changes, thus producing drastic modifications in the reflected colors. Figure 8 shows that in this case, the global appearance of the scale becomes green, as reported in [16,17]. We can also notice that the total amount of reflected energy is much lower than when the scale is in air.

 

Fig. 8. Same as Fig.7, but the scale is immersed in IPA.

Download Full Size | PDF

In the previous examples, the 16 lamellae of each elementary pattern of the tree-like structure had the same dimensions. Since it is not the case in the actual structure (Fig.2a), we investigate now how this feature affects the reflected colors. For this purpose, we will use the parameters shown in Fig.9.

 

Fig. 9. Characteristic dimensions of the second modeled scale (in nanometers). The nine lower lamellae have the same dimensions as those described in Fig.4, and the upper lamellae have decreasing widths. Apart from the bottom air layer, all the dimensions are scaled.

Download Full Size | PDF

The colors reflected by the structure of Fig.9 are represented in Figs.10 and 11. Comparison with Figs. 7 and 8 shows that the amount of reflected energy is lower than when all the lamellae have the same size. This feature appears clearly on the maps when the scale is in air. When the scale is in IPA, both left maps are so dark that the difference is not obvious, but the difference can be seen on the data from which the maps are obtained. Anyway, the overall colors are almost the same for both structures.

 

Fig. 10. Reflected colors by the structure of Fig.9 lying in air.

Download Full Size | PDF

 

Fig. 11. Reflected colors by the structure of Fig.9 lying in IPA.

Download Full Size | PDF

Another question that arises is linked to the left-right dissymmetry of the structure. Since both structures (Fig.4 or Fig.9) are not symmetric, it is interesting to know if the reflected colors are different when the incident light comes from the left or from the right. To answer this question, we have modeled a collimated beam with mean incidence (see Fig.3a) corresponding to θ ⁱ =45°, and φ ⁱ =0 (incidence from the “left”) or φ ⁱ =180° (incidence from the “right”). The incident light comes in a solid angle which spreads 10° around this mean incidence, and the reflected light is collected in a solid angle which spreads 10° around the normal (θ ^r =0). We have used the structure of Fig.4, with 16 branches for each tree, and lying in air. In this case, there is no significant difference between the reflected colors for φ ⁱ =0 and φ ⁱ =180° : they both give the same blue shade.

The reason why this structure behaves the same way when illuminated from the left or from the right is illustrated on Fig.12. Here, the structure is illuminated by a plane wave with an incidence angle θ ⁱ =45°, and we study the two cases where φ ⁱ =0 or φ ⁱ =180°. In these conditions, one reflected order escapes near the normal to the structure, and that order carries the most part of the reflected energy. The figure shows the averaged efficiencies between the two fundamental polarizations (δⁱ =0 and 90°), in order to model an unpolarized incident light. The two curves are almost identical, with a pronounced peak for the shorter wavelengths, thus generating the same blue color.

 

Fig. 12. Reflected energy near the normal when the structure with 16 branches in air is illuminated from the left (φ ⁱ =0) and from the right (φ ⁱ =180°).

Download Full Size | PDF

Finally, Fig.13 shows field maps when the structure of Fig.9 is illuminated in normal incidence by a plane wave. For λ=0.42 µm, it shows that the field does not go deeply inside the structure, and is reflected by the scale. For λ=0.60 µm, we see that the field propagates inside the structure following the dielectric regions, thus explaining a larger absorption due to the strong field value inside the lossy medium.

 

Fig. 13. Field maps of the squared modulus of the total electric field for the structure of Fig.9 illuminated by a E//polarized plane wave in normal incidence.

Download Full Size | PDF

5. Conclusion

The understanding of the iridescent colors reflected by butterflies wing is an attractive subject which has aroused several studies in the recent years. We think that the approach exposed in this paper brings two new interesting points of view.

First, the way to convert reflection coefficients into colors provides a smart technique able to bring together many numerical data into easily interpretable colormaps. Since all the content of a reflection spectrum R(λ) can be replaced by a color point in a colormap, such a colormap offers easily two free parameters in a 2D plot (and even more when using more sophisticated plottings). Such representations can be very useful to get a global insight on the reflection properties of a body.

Second, we prove that the lamellar grating theory is able to model the Morpho rhetenor scales with a great accuracy. To our knowledge, all previous studies use a simpler approach where the structure is approximated by a stack of homogeneous layers.

Using these techniques, we have recovered several properties observed experimentally, such as the color reflected by the scale lying in air or in IPA, or the spatial and spectral distributions of the energy diffracted by the wing.

At the present time, we have not considered the fact (clearly apparent in Fig. 2) that the ridges are tilted (toward the left on the Fig. 2). This asymmetry could be taken into account using our modeling, using a staircase approximation of the tilted ridges. It would lead to much larger computation times, but such a study should be affordable.

Finally, our approach could be used for various biologic species possessing an almost periodic micrometric structure.