Approches innovantes en convection thermique turbulente

Calibration des caméras : une méthode innovante

Différentes méthodes de calibration. Le premier obstacle technique au suivi de particules en trois dimensions réside dans la calibration des caméras. La projection de l’espace réel sur un capteur bidimensionnel entraine une perte d’informations qu’il faut compenser par une redondance d’informations, c’est à dire plusieurs caméras. La méthode la plus utilisée repose sur le modèle de Tsai [Tsai, 1987]. La caméra est considérée comme un capteur d’où les rayons lumineux issus du volume de mesure sont projetés par l’intermédiaire d’un diaphragme ponctuel. Bien que très efficace, cette méthode ne prend pas en compte d’éventuelles imperfections des surfaces ainsi que l’astigmatisme des dioptres plans et suppose que la lumière provient directement du milieu étudié sans aucun système optique ou dioptre avant la caméra. La calibration consiste alors à prendre une image d’une mire tridimensionnelle pour laquelle la position des points de visée est bien connue. On peut alors savoir à quel ligne dans le volume de mesure un point sur le capteur de la caméra correspond. Une méthode alternative a été développée par Svoboda et al. [Svoboda et al., 2005]. Elle consiste à suivre un point se déplaçant dans l’espace avec plus de caméras que celles utilisées pour les acquisitions. L’intérêt dans notre cas est que cela ne nécessite pas de mire à installer dans le volume de mesure qui est clos. Les caméras s’auto-calibrent alors en utilisant la redondance d’information obtenue avec le surplus de caméras. Ceci permet également d’auto-calibrer en temps réel le système optique au cours des acquisitions. Toutefois, c’est toujours le modèle de Tsaï qui est utilisé. Nous avons essayé de mettre en place cette méthode dans notre cas de figure. Elle s’est avérée inopérante. Nos caméras d’acquisition sont en effet, quelque soit leur nombre, dans un même plan de l’espace si bien que même avec des caméras supplémentaires placées hors de ce plan, l’algorithme ne permettait pas d’obtenir une calibration correcte. Une innovation : des caméras sans modèle. Une nouvelle méthode, développée au laboratoire par Miguel Lopez-Caballero, Romain Volk et Mickaël Bourgoin consiste à s’affranchir de modèle pour la caméra. La seule hypothèse est que, dans le milieu étudié, les rayons lumineux se déplacent en ligne droite 1 . Sous cette hypothèse, le trajet des rayons lumineux entre sa sortie du milieu et son arrivée sur le capteur n’a aucune importance. En théorie, si le système est stigmatique, l’image d’un point de l’espace est un point sur le capteur de la caméra. En pratique, en raison de la profondeur de champ, le lieu des points de l’espace réel où se trouve une particule détectée par la caméra constitue une ligne. Le système optique conjugue chaque pixel à une ligne de l’espace réel. La figure 3.5 illustre le propos. La méthode de calibration consiste alors à trouver, pour chaque caméra, une correspondance entre chaque pixel du capteur (qui reçoit le signal) et la ligne conjuguée dans le milieu étudié. Pour cela on utilise une mire bidimensionnelle sur laquelle sont tracés des points dont on connait précisément l’écartement. Après plusieurs essais sur des mires en papier plastifié et pour différents écartements de points, nous avons opté pour une mire usinée dans une plaque de PVC avec des points écartés de 20 mm. Cette matière assure une planéité de la mire. Les deux côtés de la plaque contiennent des points aux même positions. La figure 3.3 montre le schéma de cette mire qui est ensuite fixée sur une platine de translation Dylin SVW munie d’un compteur. L’ensemble est fixé à la cellule de convection dont on a enlevé la plaque du haut. Afin de prendre en compte les contraintes imposées aux parois en PMMA par le chauffage, la calibration est effectuée juste après les acquisitions alors que la cellule est toujours chauffée. 1. Dans le cas d’un milieu stratifié on pourrait modéliser le trajet d’un rayon lumineux s’y propageant et utiliser cette même méthode. Figure 3.3 – Plan de la plaque de calibration avec les différentes dimensions. Les trois points plus larges au centre permettent de définir l’origine et les axes qui portent X et Y . Les quatre autres à la périphérie sont des points de repères. Figure 3.4 – Images de la mire prises par deux caméras perpendiculaires. Les axes sont indiqués. Procédure de calibration. La calibration se déroule en quatre étapes. (i) Après avoir installé la mire, quelques dizaines d’images de celle-ci sont prises par chaque caméra. La mire est ensuite translatée et une nouvelle séquence d’images est prise. Cette opération est faite pour au moins trois positions. Ces positions sont centrées autour d’une valeur moyenne pour faciliter le traitement. Pour avoir la meilleure calibration possible il est préférable que la mire parcoure une large part du volume de mesure. Chaque série d’images est ensuite moyennée afin d’avoir une image par caméra et par position de la mire (voir figure 3.4). La suite du traitement est ensuite effectué sous Matlab. (ii) La seconde étape consiste, pour chaque doublet caméra/position, à calculer la transformation permettant de passer des positions en pixels sur le capteurs (Xp,Yp) aux positions dans le volume de mesure (X, Y ) ; ainsi que la transformation inverse. Pour cela, une première projection est faite en indiquant un rectangle sur les images. Cela permet une détection plus aisée des points et la correspondance de ces points à leur position réelle. Les directions −→X et −→Y sont déterminées à l’aide de trois points de repère sur la plaque permettant de fixer l’origine et les axes. Connaissant les positions dans l’espace des points présents sur le capteur, on calcule la transformation projective passant du capteur au volume de mesure. Un correction non-linéaire d’ordre trois permet de corriger, pour chaque position de la mire de calibration, d’éventuelles aberrations sphériques. (iii) La troisième partie du traitement considère chaque caméra indépendamment des autres. Les transformations effectuées précédemment permettent d’assigner à chaque pixel de la caméra une position dans le mode réel à laquelle on ajoute la position Z lue sur la platine de translation. Ainsi, sous l’hypothèse édictée plus haut, le lieu des points dans le monde réel correspondant à un même pixel pour chaque position de la mire décrit un segment de droite dont on peut déterminer un point de passage et le vecteur directeur. (iiii) Ce calcul est en fait effectué pour tout ou partie des pixels. Un interpolant est enfin construit permettant d’assigner un couple {point de passage-vecteur directeur} y compris à des valeurs de pixel non entières. L’erreur mesurée sur l’interpolant est alors de l’ordre de 0,1 mm au maximum, ce qui est inférieur au diamètre des particules que nous allons suivre. On considère ligne dans l’espace Système optique Capteur CMOS Milieu Caméra 1 Caméra 2 Particule Figure 3.5 – Schéma de principe de la correspondance entre pixel et ligne de l’espace où est situé l’objet observé. Figure 3.6 – Schéma de la correspondance tridimensionnelle d’une particule avec deux caméras. Le principe est extensible à un plus grand nombre de caméras. cela comme acceptable, d’autant plus comparé au volume de mesure. Par ailleurs, les fluctuations de température dues aux panaches entrainent des fluctuations d’indice qui viennent à l’encontre de l’hypothèse selon laquelle la lumière se déplace en ligne droite dans le milieu. Nous discuterons ce point au paragraphe 3.1.4.3. 

Suivi de particules en trois dimensions

Une fois les images acquises et les caméras calibrées, il convient de déterminer la position en trois dimensions des particules ainsi que leur trajectoire. Les routines Matlab ont été fournies gracieusement par Miguel Lopez-Caballero. Nous les avons ensuite réécrites en Python afin que le traitement puisse être parallélisé sur les fermes de calcul du Pôle Scientifique de Simulation Numérique (PSMN). On commence par effectuer une détection des particules sur les images. Elles font typiquement entre 5 et 10 pixels de diamètre. Après avoir retiré l’arrière plan, on applique aux images un filtre passe-bande afin d’éliminer le bruit et de lisser l’image, permettant ensuite une détection sous-pixel. L’image est ensuite seuillée et rendue binaire pour supprimer le bruit et détecter les centres des particules. Il faut ensuite reconstituer leur position en trois dimensions. La figure 3.7 illustre ceci. 

Correspondance à trois dimensions et suivi des particules

Nous avons vu qu’à partir des positions des particules sur le capteur, l’interpolant lié à chaque caméra nous permet de connaitre les coordonnées de la droite sur laquelle se situe la particule. La présence d’au moins deux caméras est alors indispensable. En traçant pour chaque caméra (virtuellement) la droite sur laquelle se situe la particule en question, on obtient sa position en trois dimensions par l’intersection des droites issues de chaque caméra. En réalité, on ne cherche pas une intersection mais une distance minimale entre ces droites. On fixe également un critère de distance au-delà duquel deux droites ne peuvent pas être assignées à la même particule, typiquement de 400 à 800 µm. Afin de réduire le temps de calcul, chaque droite issue de la première caméra est projetée sur les autres . Figure 3.7 – Passage de la détection sur les caméras (à gauche) à la correspondance en trois dimensions (à droite). Seule une partie des particules réellement détectées est représentée. caméras. Les particules qui se trouvent au voisinage de cette droite projetée sont sélectionnées et leur droite d’origine est tracée dans l’espace pour chercher la plus petite distance avec la droite issue de la première caméra. Ainsi on s’affranchit de la comparaison de chaque droite issue de la première caméra avec toutes celles issues des autres caméras. La figure 3.6 illustre cette technique. Toutefois celle-ci ne fonctionne pas de manière idéale. Avec seulement deux caméras et un grand nombre de particules, les situations ambigües sont nombreuses. Une droite issue d’une caméra peut passer au voisinage de plusieurs droites de l’autre caméra conduisant à des erreurs de détection. En présence de trois caméras, ces ambiguïtés sont en grande partie levées mais il en subsiste. C’est alors en essayant de suivre les particules sur deux pas de temps consécutifs que les dernières ambiguïtés sont levées. Une fois obtenues les positions en trois dimensions des particules sur l’ensemble de l’acquisition, on cherche à déterminer les trajectoires de celles-ci. Le suivi dans l’espace en trois dimensions nécessite moins de précautions que celui en deux dimensions [Ouellette et al., 2006a ; Machicoane, 2014] car les particules sont suffisamment diluées (c’est à dire suffisamment éloignées les unes des autres) par l’ajout de cette dimension supplémentaire. Aussi, au premier ordre, la recherche du plus proche voisin entre deux images consécutives est suffisante, à condition que la fréquence d’échantillonnage soit assez élevée pour que le déplacement des particules sur un pas de temps n’excède pas la moitié de la distance moyenne entre particules. En pratique, on détermine la distance maximale dmax que peut parcourir une particule dans l’intervalle de temps correspondant à cette fréquence. Pour l’obtenir on estime la vitesse maximale que peut atteindre une particule dans cette zone. Les mesures effectuées au chapitre précédent dans une expérience avec un nombre de Rayleigh plus grand indiquent que dans la zone centrale la norme de la vitesse moyenne est de l’ordre de 1 cm/s. Nous avons décidé de prendre une vitesse maximale 2 de 2 cm/s ce qui mène à dmax = 0,1 mm. On construit un ensemble de distances comprises di entre 0 mm et dmax espacées de façon logarithmique. Ensuite, pour chaque particule issue de la première image, on détermine l’ensemble des particules issues de l’image suivante qui sont dans la sphère de rayon di autour de cette particule. La particule correspondante est alors celle qui est la plus proche. Si aucune n’est trouvée, la procédure est répétée pour les particules restantes et la distance di suivante. L’utilisation de cette échelle de distances permet de limiter les erreurs de suivi qui peuvent apparaitre si l’on relie deux particules qui n’appartiennent pas à la même trajectoire. En effectuant cette procédure sur chaque paire d’images on peut reconstruire les trajectoires des particules. Toutefois, avec seulement deux caméras, cette méthode s’avère très 2. Nous verrons plus tard à l’aide des distributions de vitesses que cette valeur est adaptée. 3.1. Une expérience et des méthodes dédiées 63 limitée, et il est nécessaire d’ajouter une troisième caméra pour obtenir des trajectoires suffisamment longues.

Vitesse d’acquisition et filtrage

La convection thermique est un phénomène lent. Toutefois, il est empiriquement considéré que la résolution précise de l’accélération nécessite qu’un pas de temps entre deux images contienne au minimum 25 fois le temps de Kolmogorov [Mordant et al., 2004a]. On définit ce dernier comme : τη = rν  . (3.3) En utilisant la formule 3.2, on obtient τη = 0,39 s. La fréquence d’acquisition minimale doit alors être de 65 Hz. Toutefois, afin d’améliorer la résolution temporelle des statistiques Lagrangiennes, en particulier celle des données d’accélération, une fréquence d’acquisition de 200 Hz est préférable et a été utilisée pour les différentes acquisitions. Le sur-échantillonnage facilite de plus le filtrage et la gestion du bruit que nous allons voir. Le filtrage est en effet primordial. La vitesse et a fortiori l’accélération sont en effet des grandeurs dérivées qui amplifient le bruit sur la trajectoire. Aussi il est indispensable de filtrer les trajectoires afin de limiter l’effet du bruit sur les estimations statistiques de la vitesse et de l’accélération. Pour cela on convolue les trajectoires avec un noyau gaussien de largeur τf , chose possible car l’acquisition est sur-échantillonnée par rapport à τη. Selon que l’on veut obtenir la position, la vitesse ou l’accélération filtrée, le noyau est dérivé zéro, une ou deux fois avant la convolution. Par ailleurs, les trajectoires sont discrètes et de longueur finie. Aussi il convient de fixer une longueur au noyau gaussien. On prend ici N = 3τf . Les noyaux respectifs ont la forme suivante : — pour la position : gp(t) =    Ap exp − t 2 τ 2 f ! si |t| ≤ τf , 0 si |t| > τf , (3.4) — pour la vitesse : gv(t) =    Avt exp − t 2 τ 2 f ! si |t| ≤ τf , 0 si |t| > τf , (3.5) — pour l’accélération : ga(t) =    Aa 2 t 2 τ 2 f − 1 ! exp − t 2 τ 2 f ! + Ba si |t| ≤ τf , 0 si |t| > τf . (3.6) La constante de normalisation Ba provient de la longueur finie du filtre [Mordant et al., 2004a]. Cette dernière ainsi que Ap, Av et Aa sont déterminées à l’aide des relations suivantes (F représente le produit de convolution) :    gp(t)F1 = 1, gv(t)Ft = 1, ga(t)F1 = 0, ga(t)Ft 2 = 2. (3.7) Il convient de bien choisir la largeur du filtre pour affecter le moins possible le vrai signal. Pour cela, on s’intéresse généralement [Voth et al., 2002] à la variance de l’accélération en fonction de τf 64 Chapitre 3. Étude des statistiques Lagrangiennes par suivi de particules (figure 3.8). On observe une rapide décroissance de la variance avec la largeur du filtre avant un coude vers une décroissance bien plus lente. La première décroissance est due au filtrage du bruit tandis qu’après le coude, c’est le signal réel qui est lissé par le noyau gaussien. En extrapolant la seconde partie de la courbe par une droite, on considère que la largeur τf optimale du filtre correspond au moment où la courbe décroche de cette droite. Dans le cas présenté figure 3.8, on a donc τf = 0,29 τη.