## Abstract

We have designed a special purpose computer system for visualizing fluid flow using digital holographic particle tracking velocimetry (DHPTV). This computer contains an Field Programmble Gate Array (FPGA) chip in which a pipeline for calculating the intensity of an object from a hologram by fast Fourier transform is installed. This system can produce 100 reconstructed images from a 1024×1024-grid hologram in 3.3 sec. It is expected that this system will contribute to fluid flow analysis.

© 2008 Optical Society of America

## 1. Introduction

Analysis of fluid flow is required in various research fields, for example, aerophysics and medical engineering. Holographic particle image velocimetry systems have been studied since they have some appealing advantages. For example, they use high-quality information and allow the recording of an instantaneous three-dimensional velocity field[1, 2, 3, 4, 5].

Moreover, digital imaging technology is improving rapidly. Digital holographic techniques can easily capture the time evolution of particles using a digital camera. Members of our group have developed a complete digital holographic particle tracking velocimetry (DHPTV) system that does not use photographic films; instead it captures three-dimensional velocity vectors by using a high speed digital camera[6]. High-speed reconstruction of a particle field was successfully performed using a fast Fourier transform (FFT)[4, 5]. However, this method requires a high-speed computer to reconstruct three-dimensional images. To overcome the computational cost, we designed a special purpose computer system for DHPTV, hereafter referred to as FFT–HORN1. FFT was used for the image reconstruction. FFT–HORN1 can deal with 256×256 holograms [7]. However, in our DHPTV system, the digital camera grid size was 1024×1024. We were able to improve the FFT–HORN system to deal with more than 1024×1024 grid size holograms by utilizing an external RAM unit, FFT–HORN2. In this paper, we will report the results of the reconstruction from 1024×1024 grid size holograms.

## 2. Algorithm for DHPTV

From the Fresnel-Kirchoff diffraction equation, the equation for reconstruction is given by

where *ϕ* (*x _{i}*,

*y*,

_{i}*z*) is the amplitude of the object light at point (

_{i}*x*,

_{i}*y*,

_{i}*z*),

_{i}*λ*is the wavelength of the reference light,

*I*is the intensity of the hologram at point (

_{α}*x*,

_{α}*y*,0),

_{α}*k*is the wavenumber of the reference light, and

*N*is the number of divisions of the hologram grid in the x or y directions. Note that in Eq. (1), the hologram is in the

*z*=0 plane. Furthermore, using the Fresnel approximation under the assumption that |

*z*|≫|

_{i}*x*-

_{α}*x*|, |

_{i}*y*-

_{α}*y*|, we obtain

_{i}We introduce a new function defined by the following equation,

Substituting Eq. (4) into Eq. (3) we obtain

Equation (5) is the two-dimensional convolution integral about the *x* and *y* directions and is transformed by a two-dimensional Fourier transform,

where Φ(*n*,*m*) is the two-dimensional Fourier transform of *ϕ* (*x _{i}*,

*y*,

_{i}*z*),

_{i}*Î*(

*n*,

*m*) is the transform of

*I*, and

_{α}*G*(

*n*,

*m*) is the transform of

*g*(

*x*,

*y*).

Here, we calculate the Fourier transform of *g*(*x*,*y*) and obtain

where Δ*P* is the grid size of the hologram.

The calculation for three-dimensional image reconstruction is carried out from reproductions in several divided sections in the z-direction. Each plane is reconstructed keeping *z* constant. To calculate Eq. (3), the method described below was developed:

- calculate
*Î*(*n*,*m*) from*I*_{α} - calculate Φ(
*n*,*m*)=*Î*(*n*,*m*)×*G*(*n*,*m*) - calculate
*ϕ*(*x*,_{i}*y*,_{i}*z*) by using a two-dimensional inverse Fourier transform._{i}

From the second plane, step 1 can be omitted because *Î*(*n*,*m*) is the same when the same hologram is used.

## 3. Hardware design of FFT–HORN2

Figure 1 shows a block diagram of FFT–HORN2. This pipeline can calculate the intensity of an object using the algorithm above, with all calculations being fixed-point operations. The numbers under the slashes are the word lengths of the data over that line.

The FFT–HORN control unit controls the data flow in the pipeline to calculate the intensity of the object. The FFT CALC unit evaluates the FFT and the inverse FFT (IFFT). Here we use a Xilinx LogiCORE for the FFT CORE, which can perform a one-dimensional FFT. The G-pipeline unit evaluates Eq. (7). The Multiplier unit calculates *Î*(*n*,*m*)×*G*(*n*,*m*). The POST–Processing unit operates normalization and threshold processing. The RAM control unit controls reading and writing to the external DDR–SDRAM unit.

The hologram data (input data I in Fig. 1) are stored in the DDR–SDRAM unit through the RAM control unit. The data are sent to the FFT CALC unit for calculation of *Î*(*n*,*m*). The FFT CALC unit can perform a one-dimensional FFT. Hence one-dimensional FFT is performed twice in this pipeline to calculate the two-dimensional FFT. During the calculation, the one-dimensional FFT results are saved temporarily in the DDR-RAM unit. Next, the G-pipeline unit evaluates *G*(*n*,*m*) using the input data P1(*P*1=(λ*z _{i}*)/(

*N*

^{2}Δ

*P*)) and

*Î*(

*n*,

*m*)×

*G*(

*n*,

*m*) is calculated in the Multiplier unit. The results are sent to the FFT CALC unit and IFFT is performed in the FFT CALC unit. Lastly, the normalization of results and threshold processing are carried out in the POST–Processing unit.

We designed the FFT–HORN2 chip using a Xilinx ISE. We used the FPGA board (HORN–5 board) for FFT–HORN2. The HORN–5 board, which has been developed by Chiba university and the Institute of Physical and Chemical Research (RIKEN), has four FPGA chips, Xilinx XC2VP70–5FF1517C [9]. The chip is equivalent to about 7,000,000 gates. FFT–HORN2 connects to the host computer through a peripheral component interconnect (PCI) bus. The host computer controls FFT–HORN2, writes the intensity data of the hologram to FFT–HORN2 and sends the start signal. After FFT–HORN2 finishes calculating the intensity data of the object, the host computer reads the data.

## 4. Performance

The clock frequency of the FFT–HORN2 is 133 MHz. There is one calculation pipeline and a control unit for the calculation in a FPGA chip. We compared the calculation time of FFT–HORN2 with that of a personal computer. The specifications of the personal computer are as follows: Intel Pentium 4 3.2-GHz CPU, 2.0 GB of memory, Microsoft Windows XP operating system and Microsoft Visual C++ 2005. The FFTW subroutine library ver. 3.1.2 is used for FFT calculation [8]. Table 1 shows the calculation time for 100 reconstructed images from a hologram with a 1024×1024 grid. The calculation time of FFT–HORN2 was 3.3 sec, while the calculation time of the PC was 23.7 sec. Hence, the calculation speed of FFT–HORN2 was about 7 times faster than that of the personal computer.

Next, we show the hologram and the images reconstructed with FFT–HORN2. Figure 2 shows a 1024×1024-grid hologram of the system. Water was used for the working fluid and the injected particles were 1-micron nylon spheres. The water flow was from left to right in a micro channel of width 101.9*µm*. The grid size of the hologram (Δ*P*) was 0.4*µm*. For this simulation, holograms were generated by the system described in detail in Ref. [10].

Figure 3 shows the image reconstructed by FFT–HORN2. The reconstructions of the particles by the computer hologram algorithm were 36*µm* from the CCD surface. In Fig. 3, the black small objects are the reconstructed particles. This result shows FFT–HORN2 reconstructs the image exactly. Figure 4 is the movie of reconstructed particles flow.

Next, we reconstructed the image of a section parallel to the *x* − *z* plane (Fig. 5). Here, sections parallel to the *x*-*y* plane were reconstructed at a fixed pitch in the *z*-direction (Δ*z*=0.4*µm*) from *z*=23.2*µm* to *z*=50.8*µm* (reconstructing 70 images). The right side figure of Fig. 5 represents the *x*-*z* plane in the center of the black square region in the left side figure of Fig. 5. In this figure, the size of the particle changes in the *z* direction. This occurs because

of a special feature of the holography; namely, that the particle exists where the particle size is smallest. From this, this system can be used for the estimation of particle positions.

## 5. Conclusion and future work

We have designed and built a special purpose computer, FFT–HORN2, for DHPTV. Its clock frequency is 133 MHz and it can reconstruct 100 images from a 1,024×1,024-grid hologram in 3.3 sec. Hence, this special purpose computer system can be used for the estimation of particle positions and it is expected that this system will improve the efficiency of analysis in DHPTV.

The HORN–5 board has four FPGA chips. And the host computer which we use has four PCI bus connectors. We are devloping the FFT–HORN2 cluster system so that the intensity of an object can be calculated. In a preparatory experiment, the calculation speed of FFT–HORN2 cluster using 16 FPGA chips is about 100 times faster than that of the personal computer(Table 2). The details of this system will be described in the next our research paper. It is expected that this cluster system will result in faster reconstructions of a three-dimensional image.

## Acknowledgments

We thank Dr. T. Shimobaba and Dr. T. Sugie for their useful advice. The help of Mr. Y. Hamada is gratefully acknowledged. This research was partly supported by a Grant-in-Aid for Scientific Research (C) (20500048).

## References and links

**1. **D. H. Barnhart, R. J. Adrian, and G. C. Papen,“Phase-conjugate holographic system for high-resolution particle-image velocimetry,” Appl. Opt. **33**, 7159–7170 (1994). [CrossRef] [PubMed]

**2. **H. Memg and F. Hussain, “In-line recording and off-axis viewing technique for holographic particle velocimetry,” Appl. Opt. **34**, 1827–1840 (1995). [CrossRef]

**3. **J. Sheng, E. Malkiel, and J. Katz,“Single beam two-views holographic particle image velocimetry,” Appl. Opt. **42**, 235–250 (2003). [CrossRef] [PubMed]

**4. **U. Schnars, T. Kreis, and W. Juptner, “Direct recording of holograms by CCD target and numerical reconstruction,” Appl. Opt. **33**, 179–181 (1994). [CrossRef] [PubMed]

**5. **S. Murata and N. Yasuda, “Potential of digitalholography in particle measurement,” Opt. Laser Technol. **32**, 567–574 (2000). [CrossRef]

**6. **S. Satake, T. Kunugi, K. Sato, and T. Ito, “Digital Holographic Particle Tracking Velocimetry for 3-D Transient Flow around an Obstacle in a Narrow Channel,” Opt. Rev. **11**, 162–164 (2004).

**7. **N. Masuda, T. Ito, K. Kayama, H. Kono, S. Satake, T. Kunugi, and K. Sato, “Special purpose computer for digital holographic particle tracking velocimetry,” 14, Opt. Express , **14**, 587–592 (2006). [CrossRef]

**8. **
FFTW Home Page, http://www.fftw.org/

**9. **T. Ito, N. Masuda, K. Yoshimura, A. Shiraki, T. Shimobaba, and T. Sugie,“A special-purpose computer for electroholography HORN-5 to realize a real-time reconstruction,” Opt. Express , **13**, 1923–1932(2005). [CrossRef] [PubMed]

**10. **S. Satake, T. Kunugi, K. Sato, T. Ito, and J. Taniguchi, “Three-dimensional flow tracking in a micro channel with high time resolution using micro digital-holographic particle-tracking velocimetry,” Opt. Rev. **12**, 442–444 (2005). [CrossRef]