October 2005 

NASA/TP-2005-213931 

MOEMS Modeling Using the Geometrical 

Matrix Toolbox 

William C. Wilson 

Langley Research Center, Hampton, Virginia 

Gary M. Atkinson 

Virginia Commonwealth University, Richmond, Virginia 
The NASA STI Program Office . . . in Profile

Since its founding, NASA has been dedicated to the 

advancement of aeronautics and space science. The 

NASA Scientific and Technical Information (STI) 

Program Office plays a key part in helping NASA 

maintain this important role. 

The NASA STI Program Office is operated by 

Langley Research Center, the lead center for NASA’s 

scientific and technical information. The NASA STI 

Program Office provides access to the NASA STI 

Database, the largest collection of aeronautical and 

space science STI in the world. The Program Office is 

also NASA’s institutional mechanism for 

disseminating the results of its research and 

development activities. These results are published by 

NASA in the NASA STI Report Series, which 

includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 

completed research or a major significant phase 

of research that present the results of NASA 

programs and include extensive data or 

theoretical analysis. Includes compilations of 

significant scientific and technical data and 

information deemed to be of continuing 

reference value. NASA counterpart of peerreviewed formal professional papers, but having 

less stringent limitations on manuscript length 

and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 

and technical findings that are preliminary or of 

specialized interest, e.g., quick release reports, 

working papers, and bibliographies that contain 

minimal annotation. Does not contain extensive 

analysis. 

• CONTRACTOR REPORT. Scientific and 

technical findings by NASA-sponsored 

contractors and grantees. 

• CONFERENCE PUBLICATION. Collected 

papers from scientific and technical 

conferences, symposia, seminars, or other 

meetings sponsored or co-sponsored by NASA. 

• SPECIAL PUBLICATION. Scientific, 

technical, or historical information from NASA 

programs, projects, and missions, often 

concerned with subjects having substantial 

public interest. 

• TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and 

technical material pertinent to NASA’s mission. 

Specialized services that complement the STI 

Program Office’s diverse offerings include creating 

custom thesauri, building customized databases, 

organizing and publishing research results ... even 

providing videos. 

For more information about the NASA STI Program 

Office, see the following: 

• Access the NASA STI Program Home Page at 

http://www.sti.nasa.gov

• E-mail your question via the Internet to 

help@sti.nasa.gov 

• Fax your question to the NASA STI Help Desk 

at (301) 621-0134 

• Phone the NASA STI Help Desk at 

(301) 621-0390 

• Write to: 

 NASA STI Help Desk 

 NASA Center for AeroSpace Information 

 7121 Standard Drive 

 Hanover, MD 21076-1320
National Aeronautics and 

Space Administration 

Langley Research Center 

Hampton, Virginia 23681-2199 

 

October 2005 

NASA/TP-2005-213931 

MOEMS Modeling Using the Geometrical 

Matrix Toolbox 

William C. Wilson 

Langley Research Center, Hampton, Virginia 

Gary M. Atkinson 

Virginia Commonwealth University, Richmond, Virginia 
Available from: 

NASA Center for AeroSpace Information (CASI) National Technical Information Service (NTIS) 

7121 Standard Drive 5285 Port Royal Road 

Hanover, MD 21076-1320 Springfield, VA 22161-2171 

(301) 621-0390 (703) 605-6000 

Available electronically at http://ntrs.nasa.gov/ 

The use of trademarks or names of manufacturers in the report is for accurate reporting and does not 

constitute an official endorsement, either expressed or implied, of such products or manufacturers by the 

National Aeronautics and Space Administration. 
 1

Abstract 

New technologies such as MicroOptoElectro-Mechanical Systems 

(MOEMS) require new modeling tools. These tools must simultaneously 

model the optical, electrical, and mechanical domains and the 

interactions between these domains. To facilitate rapid prototyping of 

these new technologies an optical toolbox has been developed for 

modeling MOEMS devices. The toolbox models are constructed using 

MATLAB’s® dynamical simulator, Simulink®. Modeling toolboxes will 

allow users to focus their efforts on system design and analysis as 

opposed to developing component models. This toolbox was developed 

to facilitate rapid modeling and design of a MOEMS based laser 

ultrasonic receiver system. 

I. Introduction 

This document uses a Micro 

OptoElectroMechanical Systems (MOEMS) 

based laser ultrasonic receiver system as a 

working example to describe the modeling and 

simulation capabilities of the MOEMS toolbox. 

The system level model employs a top down 

design approach that hierarchically breaks down 

the laser ultrasonic receiver into manageable 

components. The toolbox models micro optical 

components commonly referred to as MOEMS 

devices. These component models take into 

account micro optical phenomena such as losses 

due to etch holes. Each component is modeled 

using matrix geometrical optic methods and 

image processing. The optical models employ 

geometrical system matrices for calculations 

involving both the focus and magnification factor. 

The inputs and outputs of each component consist 

of an image for dynamical manipulation and a 

state vector. The system model includes a test 

bench that simulates input stimuli and models the 

behavior of the material under test conditions. 

The components modeled include: thick lenses, 

mirrors, beam splitters, wave plates, and an 

optical detector. The wave plates are modeled 

using rotational matrices, while the other optical 

components incorporate the paraxial matrix 

methods. Other effects taken into account in the 

system model are the optical losses incurred due 

to etch-hole size and spacing. 

The primary purpose of the system model is 

not only to determine the feasibility of design but 

also to identify which component parameters are 

critical for achieving the design specifications. In 

addition, the MOEMS toolbox can also be used 

for rapid prototyping of optical systems. Example 

problems that can be quickly solved include those 

with compound lenses and multiple component 

systems. Key ideas to stress when using this 

toolbox are Geometrical Matrix Optical methods, 

top down, and hierarchical design methodologies 

and system design of MOEMS systems. 

II. Modeling Approach

The first step in top-down design methodology 

used in this example is to separate the laser 

ultrasonic MOEMS receiver system into a 

hierarchy of sub-components. The individual 

components are then modeled using a hybrid 

approach involving both geometrical matrix 

optical modeling and image processing. The 

calculations are performed by MATLAB’S®

dynamical simulator, Simulink®. Both an optical 

image and a state vector are passed to each 

component. The state vector is comprised of the 

following information: paraxial angle, paraxial 

height, phase angle, direction angle x, direction 

angle y, direction angle z, distance from 

preceding element, height of beam plane, and 

width of beam plane. Geometrical matrices are 

calculated for each component and used to modify 

the image dynamically. This information is used 

to calculate the modifications to the image and the 

new state vector that is passed on to the next 

component. 
 2

Simulink® is a dynamical simulator and as 

such can perform image processing using plane 

waves in conjunction with the geometrical matrix 

modeling. Matlab’s® Image toolbox can be used 

for manipulating an image, under the control of 

the MOEMS toolbox. Functions such as resizing 

of images, image intensity scaling, and others can 

be modified under the direction of the optical 

components discussed. 

The image processing for our initial 

application is based upon planar wave fronts, 

where each image represents the current phase of 

the wave front at that instance in time. The phase 

of the wave front is modeled using a sinusoidal 

source that modulates the intensity. This would 

allow for superposition of the waves to generate 

images such as fringe patterns. To date we have 

only used a single pixel of the image to generate 

the output signal. The detector we have modeled 

is not an array; therefore a single pixel is an 

accurate representation of the output from a 

Photo-EMF type of sensor. In the future we hope 

to create system models using the components of 

the toolbox that will perform more elaborate 

image processing. 

III. Modeling Technique 

A. Geometrical Matrices 

Geometrical matrix optical modeling 

techniques are based upon the paraxial 

assumptions, which state that the angles involved 

with ray tracing must be close to zero. In other 

words, the angles must be very close to parallel to 

the optical axis. Trigonometric functions can then 

be simplified, since sin (θ) ~ θ, and cos (θ) ~ 1, 

for very small values of θ. This approximation 

allows the tracing of rays to be accomplished 

using simple matrix manipulations. For example 

(See Fig. 1), if a ray translates a distance t1, 

through a medium with index of refraction n1, 

starting with height h, and the ray has an initial 

angle θ, then the new height h1, and new angle θ1, 

can be calculated using the following matrix 

equation. 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

=⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

h n

t

h

θ θ

1

1 0

1

1

1

1 (1) 

Fig. 1. Geometry of a ray translation. 

Similarly, the case of refraction of a light wave 

can be handled using matrices (See Fig. 2). In 

this case, if a ray translates through a medium 

with index of refraction n1, starting with height h, 

and initial angle θ, then the ray encounters a 

material with index of refraction n2 and radius of 

curvature r1, the new height h1, and new angle θ1, 

can be calculated using the following matrix 

equation. 

( )

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ − −

=⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

h r

n n

h

θ θ

0 1

1

1

2 1

1

1 (2) 

Fig. 2. Geometry of a ray experiencing 

refraction. 

The final case to investigate is when a ray is 

reflected from a surface of a material (See Fig. 3). 

In this case, if a ray is traveling through a medium 

with index of refraction n1, starting with height h, 

and initial angle θ, then the ray encounters a 

θ

Optical 

axis 

θ1 h = h1

n1 n2

Optical

axis 

h θ

θ1 h1 

t1
 3

reflective material with a radius of curvature r1, 

the new height h1, and new angle θ1, can be 

calculated using the following matrix equation. 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ − =⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

h r

n

h

θ θ

0 1

2 1

1

1

1

1 (3) 

Fig. 3. Geometry of a ray experiencing 

reflection. 

B. Component Modeling using 

Geometrical Matrices 

Using these three geometrical matrix building 

blocks, one can develop models for optical system 

components. The geometrical matrices are simply 

multiplied together to create system matrices for 

the components. Each component’s system 

matrix is used in the simulation by multiplying it 

by a subset of the current state matrix, to 

determine the new parameters for the new state 

matrix. The general equation used to calculate the 

new variables is 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅ ⎟⎟

⎠

⎞

⎜

⎜

⎝

⎛ =⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

State Vector

Current

Matrix

System

State Vector

New

_ _ . (4) 

Notice that in both figures 2 and 3 there is a 

discrepancy in the distance that the ray traverses 

to the surface of the optical element. Although 

this error is exaggerated in the figures, this error 

is assumed to be negligible by the paraxial ray 

assumptions [1]. 

1) Lenses: Because the thickness of a MOEMS 

lens is comparable to the radius of curvature of 

the lens, it is necessary to use the thick lens 

formulations and not the thin lens approximations 

[2]. See Fig. 4 for a diagram of the thick lens 

geometry. 

In the case of the thick lens the first matrix 

required for the component model is the refraction 

matrix for the incident surface of the lens. The 

second matrix is a translation matrix for the 

distance traversed within the lens. The third 

matrix describes the refraction effects as the 

optical path exits the second surface of the lens. 

Multiplying the matrices together yields the 

system matrix for a thick lens component (5). 

( ) ( )

( )

⎟

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎜

⎝

⎛

⋅

⋅ − ⋅ + ⋅

⋅ ⋅

− ⋅ + ⋅ − ⋅ + ⋅ − ⋅ ⋅⋅ + ⋅ − ⋅

=

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ −

⋅

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ −

1 2

1 2 2 2 2 1

2

2

1 2 2

2 2 2 1 2 2 2 1

2 1

2 2

2 2 2 1 2 2

1

1 2

2

2 2

2 1

0 1

1

1

1 0

0 1

1

r n

r n t n t n

n

t

r r n

r n t n t n n r

n n

r n

r n t n t n

r

n n

n

t r

n n

(5)

 Fig. 4. Thick lens Geometry 

A thin lens formulation can only be accurately 

used if the radius of curvature is much larger than 

the thickness of lens. See Fig. 5 for a diagram of 

the thin lens geometry. In the thin lens case, the 

component system matrix is reduced to simply the 

refraction matrix for the incident surface of the 

lens and the refraction matrix for the second 

surface of the lens as in (6). 

Optical 

axis

n1

r1 r2 

n2

t2

h h1 

Optical 

axis

θ

θ1 h = h1

r1
 4

Fig. 5. Thin lens Geometry 

( ) ( )

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛ − ⋅ − =

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ −

⋅

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ −

0 1

1 1 1

0 1

1

0 1

1

1 2

2 1

1

1 2

2

2 1

r r

n n

r

n n

r

n n

 (6) 

Fig. 6. Icons for the toolbox library (Lens, 

Mirror, Beam Splitter, Wave Plate, and Detector). 

The icons shown in Fig. 6 are Simulink 

hierarchical constructs that include 

implementations of the system matrices. The 

icons of the toolbox library implement basic 

optical components such as a Lens, Mirror, Beam 

Splitter, Wave Plate, and a Photo-EMF Detector. 

The Lens icon is a Simulink subsystem that is 

masked. After placing the icon, you can double 

click upon it to bring up the subsystem shown in 

Fig. 7. 

2

Lens_State

1

Thick Lens Lens_Image_Out

Matrix

Product

De m u x

2 Dem u x

State

1

Image

Paraaxial Angle

Paraaxial Height

Phase Angle

Direction Angle X

Width of the Beam Plane

Height of Beam Palane

Distance from proceeding element

Direction Angle Z

Direction Angle Y

New Paraaxial Angle

New Paraaxial Height

Fig. 7. Laser ultrasound thick lens system 

matrix. 

This subsystem illustrates the connections 

between the state vector and the image that are 

passed to the lens system matrix. Only the 

paraxial angle, paraxial height and the distance 

from the preceding element are passed to the thick 

lens system matrix. The thick lens matrix is also 

masked, and when it is clicked will bring up a 

window for typing in parameters (Fig. 8). In this 

window the lens thickness, the two radii, and the 

two indices of refraction can be input. 

Fig. 8. Thick lens parameters window. 

Optical 

axis 

n1 

r1 r2 

n2 

h h1 
 5

If you look under the mask for the thick lens 

system matrix you will see the block diagram in 

Fig. 9. In this block the distance from the 

preceding optical element (Input 2) is passed to 

the translation matrix. The translation matrix is 

multiplied by the refraction matrix for the first 

surface of the lens. The result is multiplied by the 

transfer matrix for translation that occurs between 

the lens surfaces. The resultant is multiplied by 

the refraction matrix for the second lens surface. 

The final result is the system matrix for the thick 

lens. This matrix is multiplied by the paraxial 

angle and height (Input 1) to give the new 

paraxial angle and height which are passed to the 

next optical element (Output 2). The 

magnification factor is also computed (Output 1). 

The magnification is simply the negative inverse 

of the element (2, 2), or using the ABCD matrix 

notation the magnification is -1/D. 

2

output

1

Magnification

1 0

t1 1

n1 

Transfer Matrix #1

Preceeding Optical Path

1 0

t2 1

n2 

Transfer Matrix

Internal Optical Path

-1/D

1 -(n1-n2)

 r2

0 1

Refraction Matrix

2nd Lens Surface

1 -(n2-n1)

 r1

0 1

Refraction Matrix

1st Lens Surface

X

Matrix

Multiply3

X

Matrix

Multiply2

X

Matrix

Multiply1

X

Matrix

Multiply

2

input2

1

input1

Fig. 9. Thick lens system matrix. 

In the next level up (Fig. 7) we see that the 

lens magnification which was output by the thick 

lens system matrix is multiplied by the image 

matrix to adjust the intensity. All of the 

unchanged signals, along with the new paraxial 

angle and height are passed to the next optical 

element. 

2) Mirrors: Fig. 10 shows the diagram of the 

mirror geometry. The system matrix for a mirror 

is that of a simple reflector (7) is 

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛ −

=

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎝

⎛

0 1

2 1

Matrix

System

Mirror

1

1

r

n

. (7) 

Fig. 10. Mirror Geometry 

The mirror system matrix block is masked, so 

double clicking upon the mirror icon will bring up 

the block parameter window (Fig. 11). The block 

parameter window is where mirror parameters 

such as the mirror angle, etch hole size, etch hole 

spacing, the index of refraction for the mirror 

material and the radius of curvature for the mirror 

are entered. 

Fig. 11. Mirror parameters pop up window. 

Optical 

axis

θ

θ1 h = h1

r1
 6

Underneath the subsystem mask is the mirror 

system matrix block (Fig. 12). The mirror system 

block uses the paraxial angle and height, along 

with the directional angle x. The mirror system 

matrix is multiplied by the paraxial angle and 

height to give the new paraxial angle and height. 

The direction angle x is modified to according to 

the law of reflection. The direction angle is also 

used to calculate the angle of incidence at the 

surface of the mirror which is used for de-rating 

the intensity of light reflected from the mirror 

surface. 

2

Mirror_S

1

Mirror_I_Out

Mirror

Matrix X

Matrix

Multiply

Mirror

Derating

Angle1 - u(1) + Angle1

Law of Reflection

Angle1 - u(1)

Fcn1

Dem ux

2 De m u x

Laser_S

1

Laser_I

Paraxial Angle

Paraxial height

Phase Angle

Direction angle x

Width of Beam plane

Height of Beam plane

Distance from preceding element

Direction angle z

Direction angle y

New Paraxial height

New Paraxial Angle

Fig. 12. Mirror system matrix 

The mirror matrix block (Fig. 13) determines 

the systems matrix for the mirror. In this bock the 

translation matrix for the light beam from the 

preceding optical element to the mirror is 

multiplied by the mirror reflection system matrix 

to give the overall system matrix for a mirror. 

1

Out1

1 0

t1/n1 1

Transfer Matrix

Preceding Optical Path

1 -2*n1/r1

0 1

Reflection Matrix

Mirror Surface

Matrix

Multiply

Product4

1

In1

Fig. 13. Mirror system matrix. 

Many MOEMS devices are required to have 

etch holes so that the release etchant can flow 

under the device and insure the complete removal 

of a sacrificial layer. This allows moving parts to 

be free from any substrate or unwanted 

attachments. Because etch holes have a large 

effect on the reflectivity of MOEMS mirrors they 

have been included in the models. Using data 

from [3], and lookup tables in Simulink®, a block 

has been developed that takes into account etchhole size and spacing, as well as the angle of the 

mirror when determining the intensity of light that 

is reflected from a MOEMS mirror. 

The angle of incidence at the surface of the 

mirror is passed to a mirror de-rating block (Fig. 

14). In this block the intensity of the reflected 

light is reduced using three factors. First the 

reflectance is reduced according to the angle of 

incidence. For this block look up tables are used 

to determine the outputs. The result from the 

reflectance block it multiplied by a look up table 

block that determines the reduction in light 

intensity due to the size of the etch holes in the 

mirror surface. Since almost all of the pop up 

MOEMS elements require etch holes for the 

elements to be released from the substrate, this is 

a real world issue that needs to be addressed. The 

result from this block is multiplied by the result of 

a look up table for the reduction in light intensity 

due to the spacing of the etch holes. The overall 

result of this block is multiplied by the image to 

reduce the intensity of light of the image that is 

passed to the next optical element. 

*Note only using the x angle for this calculation.

1

Derated

Reflectance vs Angle

Look-Up Table

Product2

Product1

Product

Etch Hole Spacing vs Reflectance

Look-Up Table

Spacing1

Etch Hole Spacing

Etch Hole Size vs Reflectance

Look-Up Table

Size1

Etch Hole Size

2

Angle

1

Intensity_Matrix

Fig. 14. Mirror de-rating block. 
 7

3) Wave plates: The wave plate used is a 1⁄2 λ

phase retarder and the source is linearly polarized. 

Fig. 15 shows a diagram of the wave plate 

geometry. Wave plates are used to rotate the 

phase angle of the outgoing light. 

Fig. 15. Wave Plate Geometry 

A rotation matrix like that of (8) could be used 

to determine the amount of phase rotation 

provided by a wave plate. 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

− = sin( ) cos( )

cos( ) sin( )

θ θ

θ θ

R J (8) 

However for our application which uses a half 

wave plate and assumes that the source is linearly 

polarized, the equations simplify to the level 

where matrices are not required. For this 

application we first look at the electric field as it 

is split into its x and y components as is given in 

(9). 

⎟

⎟

⎟

⎟

⎠

⎞

⎜

⎜

⎜

⎜

⎝

⎛

⋅

⋅

=⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

0

0

2

1

2

1

E

E

E

E

Y

X (9) 

Where E0 is the initial total electric field and 

EX and EY are the perpendicular and parallel 

components. These components can be used to 

calculate the angle of phase rotation for the half 

wave plate. By multiplying (9) by a matrix that 

represents linearly polarized light we get 

equations (10, 11) [4]. 

EX EX ′ = cos(θ ) ⋅ (10) 

EY EY ′ = −sin(θ )⋅ (11) 

Equations (10, 11) are used in (12) to 

determine the angle of phase rotation (α), given a 

mechanical rotation of the wave plate by angle θ. 

tan( ) cos( )

sin( ) tan( ) θ

θ

θ

α = − − = ′

′ =

X

Y

E

E (12) 

So the angle of phase rotation simplifies to 

α = −θ . (13) 

The intensity of the light (I) that has passed 

through the wave plate is found by 

2 2 ( ) ( ) EX EY I = ′ + ′ (14) 

0

2

0

2

0 sin( ) ) 2

1 cos( ) ) ( 21 I = ( ⋅ θ ⋅ E + ⋅ − θ ⋅ E = E .(15) 

The total intensity of light remains unchanged 

by a rotation of phase angle (15); therefore the 

image is passed through the phase shift block 

without any modification (Fig. 16). This block 

passes all of the signals through unchanged, 

except the phase angle. The phase angle is 

modified in the 1⁄2 wave plate sub-block. 

2

State out

1

Image Out

em

In1 Out1

1/2 Wave Plate

2

State In

1

Image In

Direction angle x

Paraxial Angle

Paraxial height

Width of Beam plane

Height of Beam plane

Distance from preceding element

Direction angle z

Direction angle y

Phase Angle

Fig. 16. Phase shifting wave plate block. 

Optical

axis

φ

h h1

φ1
 8

The half wave plate sub-block is masked, so it 

will open a block parameter window (Fig. 17). In 

this window there is only one parameter, the 

mechanical angle of rotation for the half wave 

plate. This parameter is used to adjust the relative 

intensities between the reference beam, and the 

transmission beam that has been directed at the 

surface of the material under test. 

Fig. 17. 1⁄2 Wave plate phase parameters pop 

up window. 

If you look under the mask for the half wave 

plate you will see a single function block (Fig. 

18). In this simple block the phase angle is 

multiplied by -1 to accomplish a 2*θ rotation in 

phase angle in the opposite direction. 

Laser Ultrasound

1/2 Wave Plate

1

Out1

PAngle*-1

Fcn

1

In1

Fig. 18. 1⁄2 Wave plate phase shifting block. 

4) Beam splitters: A beam splitter is modeled 

as a plane mirror with variable reflectivity and 

transmission. See Fig. 19 for a diagram of the 

beam splitter geometry. The reflectance matrix 

mentioned earlier (3) is for a curved mirror 

surface and must be modified for a plane mirror. 

To formulate the matrix for a plane mirror, you 

must take the limit of the power of the surface as 

the radius goes to infinity, which will go to zero 

[5]. 

0 ' lim 1 1 = − − →∞ r

n n

r (16) 

The zero for that element in the matrix results 

in the plane mirror reflectance matrix reducing to 

the identity matrix: 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

0 1

1 0

. (17) 

The system matrix of the component can be 

formulated (18) for the beam splitter by starting 

with a translation matrix that accounts for the 

distance and medium from the source to the plane 

mirror surface. The second matrix is the 

reflection matrix for the surface of the beam 

splitter. 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

⋅ ⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

=⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

0 1

1 0

1

1 0

1

1 0

1

1

1

1

n

t

n

t (18) 

Fig. 19. Beam Splitter Geometry 

The beam splitter sub-block (Fig. 20) has a 

single image and state vector as inputs but has 

two images and state vectors for outputs. In this 

block the phase angle is used to determine the two 

relative intensities of light that multiplied by the 

input image to create the two new images that are 

passed to the next optical elements. Also, the 

direction angle X is modified for one of the state 

vectors by adding a constant shift of 90o, while 

the other state vector direction angle X remains 

unchanged. This corresponds to the two light 

Optical 

axis

h

θ θ1

h1 t

θ2

h2 n1

φ
 9

paths that emerge from a beam splitting cube, 

where one light path passes directly through he 

cube, but the second path emerges from an 

orthogonal face of the cube. 

4

Image Out3

3

State out3

2

State out2

1

Image Out2

Product1

Product

In1

Out1

Out2

Intensity Splitter

Dem u x

pi/2

Constant

2

State In1

1

Image In1

Paraxial Angle

Paraxial height

Width of Beam plane

Height of Beam plane

Distance from preceding element

Direction angle z

Direction angle y

Phase Angle

Direction angle x

Fig. 20. Beam Splitter Block. 

The intensity splitter (Fig. 21) uses the phase 

angle to generate four different linear sections. 

The phase angle is also used to control the 

switches which select one of four linear functions 

to be output. The first function starts at zero and 

for a phase angle of zero and goes to one for a 

phase angle of 90o. The second function starts at 

one and for a phase angle of 90o and goes to zero 

for a phase angle of 180o. The third function 

starts at zero and for a phase angle of 180o and 

goes to one for a phase angle of 270o. The fourth 

function starts at one and for a phase angle of 

270o and goes to zero for a phase angle of 360o. 

2

Out2

1

Out1

Switch3

>=270o

Switch2

>=180o

Switch1

>=90o

(360-u(1))/90

Fcn4

(u(1)-180)/90

Fcn3

(180-u(1))/90

Fcn2

u(1)/90

Fcn1

1

Constant

1

In1

Fig. 21. Beam splitter intensity splitter block. 

Switch 1 outputs function 1 when the phase 

angle is less than 90o, and outputs function 2 

when the phase angle is equal or greater than 90o, 

and less than 180o. When the phase angle is less 

than 270o and is greater than or equal to 180o 

switch 2 outputs function 3. If the phase angle is 

less than 180o then the output of switch 1 is 

passed through switch 2. When the phase angle is 

360o or less and is greater than or equal to 270o 

switch 3 outputs function 4. If the phase angle is 

less than 270o then the output of switch 2 is 

passed through switch 3. In this manner all four 

functions can be output in the proper sequence, 

thus insuring that the output intensity will 

properly vary as the phase rotates through 360o. 

5) Photo-EMF Detector: Modeling of the laser 

ultrasonic receiver necessitated the modeling of 

the detector. For our application we modeled a 

Photo-EMF detector at the most rudimentary 

level. The detector block (Fig. 22) provides a 

sink for all of the state vector information in the 

form of a scope. The block receives two images 

as inputs. These are the reference beam image 

and the reflected measurement beam image. The 

reference beam image is subtracted from the 

reflected measurement beam image to give a time 

varying signal. This signal, which is roughly 

equivalent to a single pixel difference of the input 

images, is passed to the highest level of the 

Simulink model as the ultimate output from the 

system (Bottom graph of Fig. 23). 
 10

1

Analog Ouput

Scope2

Scope1

Scope

De m u x

De m u x

4

Modulated_S

3

Modulated

Intensity

2

Reference_S

1

Reference

Intensity

 

Reflected Measurement Beam

Reference Beam

Paraaxial Angle

Paraaxial Height

Phase Angle

Direction Angle X

Width of the Beam Plane

Height of Beam Palane

Distance from proceeding element

Direction Angle Z

Direction Angle Y

Paraaxial Angle

Paraaxial Height

Phase Angle

Direction Angle X

Width of the Beam Plane

Height of Beam Palane

Distance from proceeding element

Direction Angle Z

Direction Angle Y

Fig. 22. Photo-EMF detector block. 

Gaussian Constants 

The Gaussian constants can be used to 

determine a variety of properties relating to the 

component. They can be taken directly from the 

final system matrix for each component. The 

constants are found in the order shown below: 

⎟

⎟

⎠

⎞

⎜

⎜

⎝

⎛

−

−

c d

b a

 . (19) 

The principal focus point can be determined 

by: 

a

f 1 = . (20) 

The magnification factor is given by: 

0

1

b ⋅s

μ = (21) 

where 

s = −t 0 . (22) 

The image position is found to be: 

0

0

b a s

d c s S

+ ⋅

+ ⋅ = . (23) 

The principal foci are determined by: 

a

c Su

−1 = (24) 

a

b Su

− = 1

0 . (25) 

The nodal points are given by: 

a

n

n

c

S N

2

1 −

= (26) 

a

b

n

n

S N

−

= 2

1

0 . (27) 

These parameters can be used to modify the 

image (plane wave) that is being modeled in an 

application. For information on matrix optics or 

Gaussian constants refer to [1, 6]. 

C. Other Modeling Techniques 

In addition to matrix optical modeling, 

electrical, mechanical and other optical models 

can be created. New optical components such as 

micro lenses can be added to the library. Micro 

lenses from 5μm to 750μm have been fabricated 

using photoresist [7]. These lenses are roughly 

hemispherical in shape, and the surface deviation 

is on the order of one tenth of a wave [8]. For this 

application, the diameter of the photoresist lens 

can be used in our lens model as one of the lens 

radii. The index of refraction for the common 

photoresist lens material PMMA is 1.5014 [9]. 

Diffractive optical models can also be added to 

the library. In the future, mechanical and 
 11

electrical parameters will be added to the 

component models. Models can be created for 

mechanical components such as cantilever beams 

and bridges with both point and distributed loads. 

Electrical components that could be modeled 

include electrostatic gap closer actuator and 

various detectors for the system. Although we 

have not developed these components yet, they 

could easily be incorporated into the system 

framework that we have developed. 

IV. System Model of A MOEMS 

Laser Ultrasonic Receiver 

System 

A MOEMS-based laser ultrasonic receiver 

system has been chosen as the example to be used 

to investigate system modeling using the 

MOEMS toolbox. The MOEMS-based laser 

ultrasonic receiver system can be realized using 

micro optical components that are assembled on a 

silicon wafer (Fig. 24). Refer to [10, 11] for the 

basics of laser ultrasonics. This system is of 

interest for advanced nondestructive evaluation 

(NDE) technologies (for assuring the integrity of 

airframe structures, wiring, etc.) at NASA 

Langley Research Center [12, 13]. Laser 

ultrasonic inspection instruments have been 

shown to detect delaminations and fatigue cracks 

in critical structures [14]. However, laser based 

inspection instruments are generally too bulky to 

be used in the confined spaces that are typical of 

space and aircraft inspection sites. MOEMS is a 

key technology to reducing the size of current 

instruments and enabling increased inspection 

coverage in areas that were previously 

inaccessible due to instrument size and weight. 

This type of laser ultrasonic receiver is 

comprised of a laser and an optical system that 

detects the ultrasonic waves. The laser beam is 

reflected off a material in which ultrasonic waves 

have been created. The laser beam interacts with 

induced ultrasonic waves creating a reflected 

beam that is ultrasonically modulated. This 

reflected beam is directed by the micro-optical 

components to a photo-emf detector where the 

beam is combined with a reference beam to 

demodulate the ultrasonic signal [15]. The 

system is an interferometer which has a photoemf detector that outputs a signal that is 

proportional to the velocity of the fringe line 

movement. Material integrity is established by 

analysis of the demodulated signal. 

 The laser ultrasonic receiver system is 

modeled in the time domain using plane waves 

that are represented by images that are passed to 

each component. The state vector passes 

parameters to each component that are used to 

calculate any images transformations. This is 

how information is passed between the 

components of the receiver system. 

The entire laser ultrasonic receiver system is 

modeled hierarchically using Simulink®. The 

optical subsystem (Fig. 25) is one of many blocks 

and uses the MOEMS toolbox components to 

model the optical performance of the MOEMS 

devices. The optical system is divided into three 

paths; a transmission beam path on top, a 

reference beam path, and a modulated beam path 

shown on the bottom of Fig. 25. The 

transmission path has two inputs found on the 

upper left side of the diagram, the state vector, 

and an image of the initial laser beam plane wave. 

The laser beam is passed through a mirror and a 

phase shifting wave plate before entering the 

beam splitter. The reflected beam from the beam 

splitter becomes the transmission beam. The 

transmission beam path passes through a lens to 

another wave plate, to a second lens, and finally 

through a mirror and out as both an image and as 

a state vector. Both of which are sent to the 

material under test. 

The reference beam is the beam splitter’s 

transmitted beam. The reference beam is passed 

to a mirror where it is redirected to the detector. 

This beam is necessary for demodulating the 

analog signal. 

The lower set of inputs on the right side of the 

diagram represents the modulated beam path. 

This path begins with the transmission beam after 

it has reflected off of the material under test. In 

this design the beam must retrace through the 

optical path to reach the detector. At the detector 
 12

both the reference beam and the modulated beam 

are summed destructively to generate an analog 

signal. By analyzing the time of flight between 

the two pulses a determination can be made about 

the integrity of the material under test. 

Fig. 23 is a plot of the reference beam (top 

trace), the reflected measurement beam called the 

modulated beam (middle trace), and the analog 

signal generated by the detector (bottom trace). 

The sinusoidal waves represent the wave front 

phase of the plane wave images as they pass 

through a single point on the image plane. The 

red circles on the Modulated Beam graph 

highlight the areas where the ultrasonic signal is 

visibly modulated onto the carrier wave. The 

detector trace clearly shows the initial ultrasound 

pulse created within the material. It also shows 

the second pulse that is reflected from the back 

side of the material or a discontinuity within the 

material. By analyzing the time of flight between 

the two pulses a determination can be made about 

the integrity of the material under test. 

Fig. 23. Plot of the laser ultrasonic receiver 

signals 

V. Conclusion 

A geometrical matrix optical toolbox for use 

with Mathwork’s Simulink® program has been 

developed. The MOEMS toolbox contains the 

basic elements for optical system modeling; thick 

and thin lenses, a mirror, a beam splitter, a wave 

plate and a translation block. In addition to the 

optical components, it includes a block that 

models the effects of etch holes in MEMS 

mirrors. 

Modeling toolboxes allow users to spend more 

time investigating systems and design issues and 

less time developing component models. This 

toolbox was developed to facilitate rapid 

modeling of MOEMS systems; specifically a laser 

based ultrasonic receiver system being developed 

by NASA LaRC. The MOEMS toolbox can also 

be used for rapid prototyping of optical systems. 

In the future we plan on adding more 

components to the toolbox. We also intend to 

fabricate the laser ultrasonic receiver system to 

verify the simulation results generated by our 

system model. 
 13

VI. References 

[1] Gerrard, A., Burch, J. M., “Introduction to Matrix 

Methods in Optics”, John Wiley & Sons, London, 

1975, pp. 32-33. 

[2] Allen Nussbaum, Richard A. Phillips, 

“Contemporary Optics for Sci. & Eng.”, Prentice 

Hall Inc., Englewood Cliffs, 1976, pp. 5-35. 

[3] D. Goldstein, Polarized Light, 2d Ed., Chap. 11, 

(Marcel Dekker, Inc., New York, 2003, pp. 211-

240. 

[4] Warren J. Blaker, “Geometric Optics: The Matrix 

Theory”, Marcell Dekker Inc. New York, 1971 pp. 

22-24. 

[5] J. Alda, "Paraxial Optics", Encyclopedia of 

Optical Engineering", Driggers, R., Ed., Marcel 

Dekker Inc., New York, 2003, pp. 1920-1931. 

[6] Jun Zou; M. Balberg, C. Byrne, Liu Chang; D.J. 

Brady, Optical properties of surface 

micromachined mirrors with etch holes, MEMS, 

Journal of , Vol. 8, Issue. 4 , Dec. 1999, pp. 506 -

513 

[7] D. Daly, R. Stevens, M. Hutley, N. Davies, “The 

Manufacture of Microlenses by Melting 

Photoresist”, MicroLens Arrays, IOP Pub. Ltd., 

Bristol, pp. 23-34. 

[8] T. R. Jay, M. B. Stern, and R. E. Knowlden, 

‘‘Effect of refractive microlens array fabrication 

parameters on optical quality,’’ SPIE Proc. Vol. 

1751, Miniature & Micro-Optics: Fabrication and 

System Applications II, July 19-20 1992, San 

Diego, pp. 236-245 

[9] Gregory T. A. Kovacs, “Micromachined 

Transducers Sourcebook”, New York: McGrawHill, 1998, pp. 490. 

[10]Jean-Pierre Monchalin, Christian Néron, Jean F. 

Bussière "Laser-ultrasonics: from the laboratory to 

the shop floor", Physics in Canada, March, Vol. 51 

No 2, http://www.ultraoptec.com/press/j-pph.html. 

[11]M. Klein, B. Pouet, P. Mitchell, “Photo-emf 

detector enables laser ultrasonic receiver”, 

Detectors, Laser Focus World, Aug. 2000, url: 

http://www.newport.com/Support/Magazine 

Features/00390-687.pdf 

[12]Charles E. Harris, James H. Starnes, Jr., and Mark 

J. Shuart, “An Assessment of the State-of-the-Art 

in the Design and Manufacturing of Large 

Composite Structures for Aerospace Vehicles”, 

NASA Technical Memorandum NASA/TM-2001-

210844, Langley Research Center, Hampton, 

Virginia, April 2001, pages 50. 

[13]W. H. Prosser, T. L. Brown, S. E. Woodard, G. A. 

Fleming, and E. G. Cooper, “Sensor Technology 

for Integrated Vehicle Health Management of 

Aerospace Vehicles”, AIP Conference 

Proceedings, March 27, 2003, Volume 657, Issue 

1, pp. 1582-1589. 

[14]C. B. Scruby, L. E. Drain, “Laser Ultrasonics: 

Techniques and Applications”, Adam Hilgar, 

Bristol, 1990, pp. 325-353. 

[15]David M. Pepper, “The Emergence of Nonlinear 

Optics in the Factory”, Optics and Photonics News 

Magazine, May 1997. 
 14

Fig.24. A MOEMS Laser Ultrasonic Receiver System 

Silicon Substrate 

Detector

Mirror 

M1 

Mirror 

M3 

Beam 

Splitter 

BS1

Polarizing 

Grating 

P2

1⁄2 

λ

Wave 

Plate P1 

θ

Lens 

L1 

Lens 

L2 

Mirror 

M2 

Material Under Test 

Transmission Beam PathModulated Beam PathReference Beam Path
 15

Fig.25. Block Diagram of Optical Subsystem 

Laser

Analog Output Signal

Modulated Beam Path

Transmission Beam Path

Reference Beam Path
REPORT DOCUMENTATION PAGE Form Approved

OMB No. 0704-0188

2. REPORT TYPE 

Technical Publication

 4. TITLE AND SUBTITLE

MOEMS Modeling Using the Geometrical Matrix Toolbox

5a. CONTRACT NUMBER

 6. AUTHOR(S)

Wilson, William C.; and Atkinson, Gary M.

 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

NASA Langley Research Center

Hampton, VA 23681-2199

 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

Washington, DC 20546-0001

 8. PERFORMING ORGANIZATION

 REPORT NUMBER

L-19179

10. SPONSOR/MONITOR'S ACRONYM(S)

NASA

13. SUPPLEMENTARY NOTES

An electronic version can be found at http://ntrs.nasa.gov

12. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified - Unlimited

Subject Category 33

Availability: NASA CASI (301) 621-0390

19a. NAME OF RESPONSIBLE PERSON

STI Help Desk (email: help@sti.nasa.gov)

14. ABSTRACT

New technologies such as MicroOptoElectroMechanical Systems (MOEMS) require new modeling tools. These tools must

simultaneously model the optical, electrical, and mechanical domains and the interactions between these domains. To facilitate

rapid prototyping of these new technologies an optical toolbox has been developed for modeling MOEMS devices. The

toolbox models are constructed using MATLAB's dynamical simulator, Simulink. Modeling toolboxes will allow users to

focus their efforts on system design and analysis as opposed to developing component models. This toolbox was developed to

facilitate rapid modeling and design of a MOEMS based laser ultrasonic receiver system.

15. SUBJECT TERMS

EMS; Geometrical Optics; MOEMS; Optical Modeling; Ultrasonic Transducer

18. NUMBER

 OF 

 PAGES

20

19b. TELEPHONE NUMBER (Include area code)

(301) 621-0390

a. REPORT

U

c. THIS PAGE

U

b. ABSTRACT

U

17. LIMITATION OF 

 ABSTRACT

UU

Prescribed by ANSI Std. Z39.18

Standard Form 298 (Rev. 8-98)

3. DATES COVERED (From - To)

5b. GRANT NUMBER

5c. PROGRAM ELEMENT NUMBER

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

23-064-30-12

11. SPONSOR/MONITOR'S REPORT

 NUMBER(S)

NASA/TP-2005-213931

16. SECURITY CLASSIFICATION OF:

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, 

gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this 

collection of information, including suggestions for reducing this burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and 

Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person 

shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number.

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS.

1. REPORT DATE (DD-MM-YYYY)

01- 10 - 2005