Table Of ContentMicron-scale Light Transport Decomposition Using Interferometry:
Supplementary Material
IoannisGkioulekas AnatLevin Fre´doDurand ToddZickler
HarvardUniversity WeizmannInstitute MassachusettsInstituteofTechnology HarvardUniversity
Abstract image-basedrendering,imageediting,andmeasuringsceneshape
inthepresenceoftranslucencyandinterreflections.
Wepresentacomputationalimagingsystem,inspiredbytheopti-
EachelementofthetransportmatrixinEquation(1)stillrepresents
calcoherencetomography(OCT)framework,thatusesinterferom-
asumofdifferentscenepaths.AsdepictedbythebluepathsofFig-
etry to produce decompositions of light transport in small scenes
ure1,sub-surfacescatteringandinterreflectionsmeanthatthereare
orvolumes.Thesystemdecomposestransportaccordingtovarious
typicallymanypathsthatoriginateatthesamesourceelementand
attributesofthepathsthatphotonstravelthroughthescene,includ-
arriveatthesamecameraelementbuttakedifferentroutesthrough
ingwhereonthesourcethepathsoriginate,theirpathlengthsfrom
thescene. Computationalimagingcanbeusedtoanalyzetheseas
sourcetocamerathroughthescene,theirwavelength,andtheirpo-
well, by decomposing each entry T of the transport matrix ac-
larization. Sinceitusesinterference,thesystemcanachievehigh pl
cordingtothecontributionsfromphotonpathsthathavedifferent
pathlengthresolutions,withtheabilitytodistinguishpathswhose
opticallengthsτ. Thisleadstothenotionofapathlength-resolved
lengthsdifferbyaslittleastenmicrons. Wedescribehowtocon-
lighttransportmatrixTτ, τ ∈{τ ,...,τ },wherepathlengths
structandoptimizeanopticalassemblyforthistechnique,andwe min max
arediscretizedintoafinitesetof∆τ-sizedlengthintervals.Thefull
buildaprototypetomeasureandvisualizethree-dimensionalshape,
lighttransportmatrixcanthenbewrittenas
directandindirectreflectioncomponents,andpropertiesofscatter-
ing,refractive/dispersive,andbirefringentmaterials.
T=(cid:88)Tτ. (2)
CRCategories: I.4.1[ImageProcessingandComputerVision]: τ
DigitizationandImageCapture—Imaginggeometry;
EachentryTτ isthefractionofphotonsoriginatingatthelthsource
pl
elementthatarriveatthepthcameraelementhavingtraveledpaths
Keywords: lighttransport,interference,waveoptics
whose optical lengths are in the interval τ ±∆τ/2. There is a
growingnumberofmethodsformeasuringeitherentireprojections
1 Introduction
orisolatedelementsofthepathlength-resolvedtransportmatrix.In
addition to enhancing the applications listed above, these decom-
Inimaging,photonsleaveasource,travelthroughascene,andare
positionscanbeusedtovisualizelight-in-flightthroughtable-top
collectedbyacamera. Aconventionalimagemeasuresthesumof
scenes,andfor“imagingaroundcorners.”
allphotonsthatarriveateachcamerapixel,regardlessofwhereon
the source they originate, or which path they travel from there to Weintroduceanewcomputationalimagingsystemthatusesopti-
thecamera(Figure1).Thissummationconflatesinformationabout cal interferometry to produce high-fidelity light transport decom-
theshapesandmaterialsinthescene. Inthispaper,weshowhow positions. Our system uses optical configurations that are varia-
interferometrycanbeusedtomeasuredecompositionsoftheseper- tions of the classical Michelson interferometer, and our analysis
pixelsums,distinguishingbetweenphotonpathsthatdifferintheir buildsontechniquesthathavebeenusedforopticalcoherenceto-
lengthandlocationoforigin. mography(OCT).Oursystemiscomplementarytoexistingcompu-
tationalphotographymethodsforproducingsuchdecompositions,
In discrete terms, where the camera sensor is partitioned into P
excellinginusecaseswhereitisnecessarytoimagesmallscenes
area elements (“pixels”) and the source is partitioned into L area
at very high spatial and pathlength resolutions. In particular, the
elements, the energy that is transferred from source to camera is
useofinterferometryallowsoursystemtoachievepathlengthres-
oftendescribedbyaP×LmatrixTcalledthelighttransportma-
olutionsaslowas10µm,whichisnecessarytoanalyzetransport
trix.EachentryT inthismatrixrepresentsthefractionofphotons
pl eventscausedbymaterialeffectslikedispersionandscattering.
followingpathsoriginatingatthelthsourceelementandarrivingat
thepthcameraelement. WithreferencetoFigure1,orangeversus OurpaperbeginswithbackgroundontheMichelsoninterferome-
bluepathswouldcorrespondtoentriesTpl fordifferentvaluesof terandthenotionsofspatialandtemporalcoherencelength. We
sourcelocationl. Conventionalimagingsimplymeasurestheim- then introduce a mathematical model of interferometry in terms
agei(aP-vector)forsomepatternl(anL-vector)projectedfrom of(acontinuousversionof)thepathlength-resolvedlighttransport
thesourceaccordingtothelighttransportequation[Ngetal.2003] matrix. We use it to show how sources with different coherence
propertiesenabledifferentkindsoflighttransportdecompositions,
i=Tl. (1) differentiatinglightpathsintermsoftheirendpointlocations,their
opticallengths, orcombinationsofthesetwo. Wealsocharacter-
There are a variety of computational imaging methods for mea-
izeresolutionandnoiseperformance, andpresentaperformance-
suring elements of the transport matrix. This has applications in
optimized optical design that additionally allows resolving trans-
port in terms of wavelength and polarization. Our prototype has
threespectralchannels,twopolarizationchannels,aworkingvol-
ume of 2cm H×2cm W×1cm D, and spatial and pathlength
resolutions that are both 10µm. We use this prototype to obtain
micron-scale decompositions of light transport in scenes contain-
ingreflection,refraction,dispersion,scattering,andbirefringence.
2 Related Work
Lighttransportdecomposition. Methodsfordecomposinglight
transport, both spatially and in terms of pathlength, can be de-
scribedascapturingdifferentimagesiofastaticsceneunderuni-
adaptive function w(τ) that selects for each such non-zero entry
mirror theshortestpathlengthτ forwhichtransportisnon-zero. Inthese
reflection cases, thepathlengthresolution∆τ determineshowcleanthere-
sulting separations are, with sharper resolutions allowing less of
theglobalcomponent—suchasthatcausedbysub-surfacescatter-
inginFigure1—tobleedintothedirectone.Secondisthetemporal
direct
paths diffusewall frequencyprobingmethodofO’Tooleetal.[2014b]. Thisallows
decompositionswithmoregeneralchoicesofmatrixMandfunc-
tionw(τ),forinstancecapturingtransientimagescorrespondingto
fixedspatialprobingpatternsM.
backscattering
Interferometry. Optical interferometry refers to a large set of
Figure 1: Differentiating light paths. A pixel in a conventional imagingtechniquesthatexploitinterferencebetweenoneormore
image measures the sum of all lights paths arriving at that pixel. electromagneticfields[Hariharan2003].Despitetheirpopularityin
Entriesofalighttransportmatrixdistinguishpathswithdifferent otherdisciplines,interferometrictechniqueshaverarelybeenused
startinglocations(orangevs.blue),andentriesofa“pathlength- forcomputergraphics. AnexceptionisCossairtetal.’s[2014]use
resolved” light transport matrix additionally distinguish paths of ofaMichelsoninterferometerforrefocusing.
differentlengths(e.g.,lightvs.darkblue).
Many interferometric techniques can be interpreted using Equa-
formilluminationl=1accordingtotheequation tion (3). Abramson [1983] used holographic techniques to cre-
atevisualizationsoflight-in-flightsimilartorecenttransientimag-
i=(cid:88)w(τ)(M(cid:12)Tτ)1. (3) ing; andtheconnectionbetweeninterferometryandtime-of-flight
sensorshasbeendiscussedelsewhere[Schwarteetal.1995;Luan
τ 2001]. Here,were-establishthisrelationship,andweexpanditto
relateinterferometrytoothertypesoflightpathdecompositions.
Thisisapathlength-resolvedvariantofthetransportprobingequa-
tionofO’Tooleetal.[2012]. Theoperator(cid:12)representspointwise
Amonginterferometrictechniques,theonemostcloselyrelatedis
multiplication(Hadamardproduct). MisaP ×Lbinarymatrix:
optical coherence tomography (OCT) [Huang et al. 1991], which
setting M to zero removes contributions of paths beginning at
pl is used for range scanning and visualizing volume cross-sections
pointl onthesourceandarrivingatpixelponthesensor. Simi-
atmicronresolutions. Therearemanyvariantsbasedonsingleor
larly, w(τ)isascalarbinaryfunctionthatcanbeusedtoremove
multipleshots(Fourier-domainvs.time-domainOCT)andsingle-
contributionsfrompathsoflengthτ. Inconventionalimaging,M
orfull-fieldillumination. Inourpaper,weadaptafull-field,time-
andw(τ)areequalto1everywhere, andEquation(3)reducesto
domainOCTconfigurationforbroadercomputationalimaging.We
thestandardlighttransportEquation(1)withuniformillumination.
re-interpret conventional OCT measurements using Equation (3),
Computationalimagingmethodsforlighttransportdecomposition
andthenbuildonthistomeasureotherlighttransportdecomposi-
collectmeasurementswithdifferentMandw(τ).WewilluseFig-
tionswithdifferentMandw(τ)choices.Wefocusontime-domain
ure1tovisualizepathscapturedbyvariousapproaches.
OCT,asopposedtoFourier-domain,becauseitadditionallyallows
SpatialdecompositionmethodsinducedifferentmatricesMwhile decomposingtransportaccordingtowavelengthandpolarization.
keepingw(τ) = 1everywhere. WithreferencetoFigure1,these
Comparison and trade-offs. The main advantage of using in-
methodscanbeusedtorecordcontributionsfromonlytheblueor
terferometrytodecomposetransportishighpathlengthresolution.
onlytheorangepaths. Byusingdifferentarrangementsofacam-
Fempto-secondlasersandtime-of-flightsensorscanachievepath-
era and a projector, and different sensor and illumination coding
lengthresolutionsofabout600µmand10,000µm, respectively.
strategies,thesemethodscaptureindividualelementsofthetrans-
Bycomparison,weachievepathlengthresolutionsofabout10µm.
portmatrix,orsparseorlow-rankapproximationstoit[Senetal.
Otherinterferometrictechniquescouldbeusedtoproduceresolu-
2005;Peersetal.2009;Wangetal.2009;Baietal.2010;O’Toole
tionsthatareevenhigher,reachingsub-micronscales,forinstance
andKutulakos2010;O’Tooleetal.2012;O’Tooleetal.2014a].
usingverywide-bandsourcesandpost-captureprocessingthatin-
There are also methods that decompose according to pathlength cludes phase shifting. However, this would come at the expense
only,usingamatrixM = 1andcontrollablefunctionsw(τ). Re- of the ability to probe different spectral channels. Since we are
ferringtoFigure1, thisallowsmeasuringthesuperpositionofall interestedinmeasuringwavelength-dependenteffects,suchasdis-
paths(blueandorange)thathavethesamelength. Thishasbeen persion, birefringence, and scattering, we do not attempt to push
dubbedtransientimagingandcanbeusedtoproducelight-in-flight pathlengthresolutiontothoselimitshere.
visualizationsofphotonspropagatinginascene. Implementations
include using a combination of pulsed laser and ultra-fast detec- This increased pathlength resolution is valuable in cases where it
tor [Velten et al. 2013; Wu et al. 2014b; Wu et al. 2014a], or a isnecessarytoobtainhigh-fidelityscansofsmallvolumes,suchas
combinationofsourceandcamera(time-of-flightsensor)thatare inlaboratorymaterialmeasurementsorinsceneswithnaturalmi-
synchronouslymodulatedatradio-frequenciesintime[Heideetal. crostructures. Whilemacrophotographyreadilyprovidesawayto
2013;Kadambietal.2013;O’Tooleetal.2014b;Heideetal.2014]. imagesuchscenesatspatialresolutionsthatapproachthediffrac-
Inthelattercase,theperformancecriticallydependsonthemethod tionlimit(afewmicrons),comparablepathlengthresolutionscan-
usedtosolvetheso-called“multi-pathinterference”problem[Dor- notbeachievedwithouttheuseofinterferometry.
rington et al. 2011; Heide et al. 2014]. The equivalence of path-
The main disadvantages of using interferometry to decompose
lengthdecompositionandlightpropagationhasalsobeenexploited
transportarereductionsinfieldofview,depthoffield,andspeedof
forefficientrenderingoflight-in-flight[Jaraboetal.2014].
capture.Ourworkingvolumecanonlyaccommodatematerialsam-
Finally, there exist methods that provide simultaneous decompo- plesorsmallobjects,andoursystemisnotappropriatefortable-top
sitions of space and pathlength. In the computer graphics liter- orroom-scalescenes. AsweexplaininSection3,ourmethodop-
ature, these fall into two categories. First are methods that de- erates by slicing a volume in micron-sized depth intervals. This
compose transport into a direct component (“one-bounce paths”) meansthatmeasuringa1cmdeepvolumerequires10,000images
and an indirect or global component (“multi-bounce paths”) [Na- andseveralhours,duringwhichthescenemustremainstill. This
yaretal.2006;Gupta,M.andAgrawal,A.andVeeraraghavan,A. isexacerbatedbythefactthatouropticalsetupisverysensitiveto
and Narasimhan, S.G. 2011; Reddy et al. 2012]. These imply a vibrationsinducedbyenvironmentsources,andevenmicron-scale
matrixMthathasonlyonenon-zeroentryperrow,andaspatially- vibrationscanseverelyaffectourmeasurements.
Ms scene scene z x1 oγo1γ2 ooγ4γ3 x1 τ
d Mr Mr zo uin,θ us,θ uin,θ ur,τ,θ x
source s d(cid:48) d
r r x
(a) (b) (c) (d) (e)
iii.obliqueillumination
beamsplitter Figure3:NotationandcoordinatesystemusedinSections3and4.
(a)ThepartsoftheMichelsoninterferometercorrespondingtothe
camera
scene and reference arms. (b-e) The scene and reference arms
shown unfolded in the same coordinate system. (b) and (c) show
i:d(cid:48) ii:d
r r thetargetarms,withlightingandcamerarespectively. (d)and(e)
Figure 2: Michelson interferometer. An input beam is split by a
showthereferencearms,withlightingandcamerarespectively.
beamsplitterintotwocopiesthattraversetomirrorsM andM at
r s
distancesd andd fromthesplitter,respectively. Thetwocopies
r s M .Bydetectinginterferencepatternsinthemeasuredimages,we
reflectbackandrecombineatthebeamsplitterbeforebeingimaged. r
canseparatelightpathsofdifferentpathlengths,uptoaresolution
Insets(i)and(ii)showimagesfortwopositionsofmirrorMr that equal to the temporal coherence length. We can additionally use
eitherinduceinterference(d(cid:48)r ≈ds)ordonot(dr (cid:54)=ds).Inset(iii) sourceswithsmallspatialcoherencetoseparatepathsbylocation.
showsasetupthatcanbeusedtoachieveobliqueillumination. Thisisbecauseonlylight-pathsthatbeginsufficientlyclosetoeach
otheronthesource(i.e.,arelessthanonespatialcoherencelength
Anotherlimitationofourapproachrelativetonon-interferometric apart) and combine at the same sensor pixel will interfere. Note
techniquesisreducedflexibilityinthespatialprobingpatternsM. that this discussion is based on the coaxial interferometer design
Wecontrolthesepatternsbyusingmirrorsofdifferentshapes(Fig- ofFigure2, whichisthedesignweuseexclusivelyinthispaper.
ure 6). This means that we can only induce patterns that corre- Obliqueilluminationisachievable,butitwouldrequireadditional
spondtomanufacturablemirrors,andthatchangingthepatternM componentsasshowninFigure2.iii.
requires substantial manual reconfiguration. This is in stark con-
trasttospatialprobingmethods[O’Tooleetal.2012;O’Tooleetal. Thisdiscussionhighlightstheimportanceoftemporalandspatial
2014a;O’Tooleetal.2014b]thatallowcomplexpatternsandsim- coherence,whichwewillnowdescribeanalytically.Theyarewave
plereconfigurationinsoftware. phenomena, butscalarwavetheoryissufficienttodescribethem,
without having to employ the full complexity of electromagnetic
The third limitation of our approach is shared by all methods for theory. To simplify notation, and without loss of generality, we
pathlengthdecomposition,namely,thereducedsignal-to-noisera- consideratwo-dimensionalworld,withspatialpointsrepresented
tiowhenmultiplepathswithwidelyvaryingenergiesarriveatthe bytheir(x,z)coordinates.Figure3showsthecoordinatesystem.
samepixel. Inextremecases,thismakesitdifficultorimpossible
toresolvelow-energypaths,andunlikeconventionalHDRphotog- 3.1 TemporalCoherence
raphy,thiscannotbesolvedbycapturingmultipleexposures.
Webeginbyassumingthatthesourceproducesaplanewaveprop-
3 Optics Background agating in the z direction, so that the electromagnetic field is in-
dependentofx. Suchawavecanbeproducedbytransferringthe
Wefirstpresentthebackgroundnecessaryforunderstandinginter- outputofapointsourcethroughalens,asinFigure8. Ifthewave
ferometry and our interpretation of it as decomposing transport. ismonochromaticatwavelengthλ,theelectromagneticfieldatpo-
Our description of temporal coherence follows Goodman [2000], sition(x,z)andtimetcanbedescribedas
andourdescriptionofspatialcoherencefollowsLevinetal.[2013].
u(x,z,t)=a·exp(ik(−zη+ct)), (4)
InterferometrybeginswiththeMichelsoninterferometer,shownin
Figure2,withalightsourceofsimultaneouslysmalltemporalco- wherecisthespeedoflightinvacuum,ηistherefractiveindexof
herenceandspatialcoherence. Wewillexplainthesenotionsand themedium,k = 2π/λisthewavenumber,andaisthecomplex
theoperationoftheMichelsoninterferometerintuitively,andthen amplitudeofthesource.
describe them analytically. The setup uses a beamsplitter to split
thelightwaveemittedfromthelightsourceintwoparts. Onepart Real waves are not perfectly monochromatic, but are better de-
istransmittedtowardsareferencemirrorMrplacedonatranslation scribedassuperpositionsofplanewavesofdifferentwavelengths,
stage.Theotherpartisreflectedtowardsthetargetscene,whichfor
(cid:90)
ourintuitiveexplanationisanothermirrorM ,butcangenerallybe
s u(x,z,t)= a ·exp(ik(−zη+ct)) dk, (5)
anyarbitraryscene.Werefertothesidesofthesetupcontainingthe k
k
referenceandtargetmirrorasthereferencearmandtargetarm,re-
spectively.Bothwavesarethenreflected,recombinedatthebeam- where the complex number ak is the amplitude of each compo-
splitter,andimagedbyacamera.Wedenotebyd ,d thedistance nent. Typicalrealworldsourcescanbemodeled[Goodman2000]
r s
betweenthebeamsplitterandmirrorsMrandMs,respectively. bfoyrmaslsyuamnidngpothwaetrask|aare|2rsaanmdopmledvafrroiamblaesGwaiutshspiahnasweisthsammepalnedk¯uannid-
k
Thesmalltemporalcoherenceofthelightsourcemeansthat,when standarddeviation∆ ,
k
|d −d | is larger than an amount called the temporal coherence
r s
leernagwthilolfmtehaesusoreuracne,imthaegteweoqwuaalvteostdhoensuomtinotfertfheeretwanodimthaegceasmo-f |ak|2 ∝e−(k2−∆k¯2k)2. (6)
the mirrors, as shown in Figure 2.ii. If instead we use the trans-
lationstagetorepositionMr toapositiond(cid:48)r sothat|ds−d(cid:48)r|is Thestandarddeviation∆kisthespectralbandwidthofthesource:
smallerthanthecoherencelength,thetwowaveswillinterfereand smallervaluesof∆kindicateamoremonochromaticsource.
thecamerawillmeasureafringepattern,asshowninFigure2.i.
NowsupposeasensorrecordsatpixelxameasurementI(x)ofthe
WhenwereplacemirrorM withageneralscene,wecancapturea superpositionoftwocopiesofthewaveufromEquation(5)that
s
sequenceofimagesatdifferenttranslationsofthereferencemirror havetraveleddifferentdistanceszandz+τ. Themeasurementis
averagedoverexposuretime,producing thephaseshift,resultinginEquation(10).(cid:3)
I(x)=(cid:10)|u(x,z,t)+u(x,z+τ,t)|2(cid:11)t (7) 3.2 SpatialCoherence
=2I (x)+2Re{corr(x,τ)}, (8)
o
Sofarwehaveconsideredthecaseofaperfectplanewavepropa-
where(cid:104)·(cid:105) denotesexpectationovertime,I (x)=(cid:10)|u(x,z,t)|2(cid:11) gatingalongthezaxis.Suchawavehasperfectspatialcorrelation,
istheintetnsityoftheoriginalwave,and o t inthesensethatifweknowthefieldvalueu(x,z,t),wecanalso
predictitsvalueatanyspatialshiftu(x+ξ,z,t). Thisisonlyre-
corr(x,τ)=(cid:104)u(x,z,t)∗·u(x,z+τ,t)(cid:105) (9) alizablebyusingalenstocollimatetheoutputofasourcethathas
t infinitesimalarea. However, manyrealsourceshavelargereffec-
is the correlation of the two waves. If zero correlation exists be- tive areas. This induces another effect known as spatial incoher-
tweenuanditsshiftedcopy,themeasurementwillsimplybeequal ence[Levinetal.2013],whichalsoaffectsinterferencecontrast.
to twice the intensity of the original wave, 2I (x). If positive
o
correlation (constructive interference) or negative correlation (de- Eachpointonanareasourceemitsanindependentwave,resulting
structive interference) exists, the measurement will vary between cumulativelyinacollectionofindependentplanewavesu (x,z,t)
θ
0and4Io(x),dependingontherelativeshiftτ ofthetwocopies. overasmallangularrangeθ ∈ [−Θ/2,Θ/2]. Asintheprevious
Themaximumvalueisachievedwhenthereisperfectcorrelation section,weconsidertwocopiesofthecumulativewave,butinthis
(τ =0).Wewillrefertotherealpartofthefieldcorrelation,equal casewithonecopyshiftednotonlyinthez,butalsointhexdirec-
tothedifferenceI(x)−2Io(x),astheinterferencecontrast. The tion. Theintensityatapixelonthesensorisgivenbyintegrating
temporalcoherencelengthofthesource,whichwedenotebyLc, theintensitiesfromEquation(7)forallindependentwaves,
isthelargestvalueofτ forwhichsignificantcorrelationexistsand
thereforesignificantinterferencecontrastismeasured. Toevaluate (cid:90) Θ/2
thetemporalcoherencelength,weusethefollowingclaim. I(x)= (cid:10)|uθ(x,z,t)+uθ(x+ξ,z+τ,t)|2(cid:11)t dθ. (14)
−Θ/2
Claim1. ThefieldcorrelationmagnitudeisaGaussianwithstan-
darddeviationinverselyproportionaltothespectralbandwidth∆ , Thisimpliesthatthecorrelationis,
k
corr(x,τ)=eik¯τG∆−k1(τ), (10) corr(x,ξ,τ)=(cid:90) Θ/2(cid:104)uθ(x,z,t)∗·uθ(x+ξ,z+τ,t)(cid:105)tdθ, (15)
−Θ/2
where
G∆−1(τ)=e−21(∆kτ)2. (11) andsimilarlyforthefieldintensityIo(x)andtheinterferencecon-
k trastI(x)−2I (x). Toevaluatehowfasttheinterferencecontrast
o
decaysasafunctionofspatialshift,weusethefollowingclaim.
Beforeprovingtheaboveproperty, weconsiderandexperimental
visualization of the result. Figure 4 shows the interference com- Claim 2. For an area light source with an angular range θ ∈
ponent we measured using the Michelson interferometer setup of [−Θ/2,Θ/2],thecorrelationdecaysas
Figure2,plottedasafunctionofpathlengthdifference. Asshown
in the close-up plot, the profile looks like a fringe pattern, which corr(x,ξ,τ)≈W (ξ)eik¯(τ)G (τ), (16)
variesasafunctionofthephasedifferenceintroducedbyτ. When ∆c ∆−1
k
k¯τ = (2p + 1)π (p any integer), we observe destructive inter-
ference(zerointensity); andwhenk¯τ = (2p)π, weobservecon- where (cid:18) ξ (cid:19) λ¯
structiveinterference(twicethesignalintensity). Asshowninthe W∆c(ξ)=sinc ∆ , ∆c = 2Θ. (17)
zoomed-outplot,thecontrastofthefringepatterndecaysandeven- c
tually becomes almost zero once the pathlength difference τ ex-
ceeds some amount. We can define this amount as the temporal
Weseethat,forareasources,theinterferencecontrastdecaysvery
coherencelengthofthesource,whichfromClaim1isthestandard fast as the spatial shift ξ increases, and no significant contrast is
deviationofthecorrelationwindow,Lc =1/∆k. Weobservethat measured for values of ξ larger than the quantity ∆ of Equa-
Lc isinverselyproportionaltospectralbandwidth. Thelargerthe tion(17).Analogouslytotheprevioussection,wedefince∆ asthe
c
bandwidthofthesource,thesmallerthetemporalcoherencelength, spatialcoherencelength. Thespatialcoherencelengthisinversely
andthereforethehighertheresolutionatwhichwecandiscriminate
proportional to the angular extent of the source. The wider the
betweenpathsofdifferentlengths.
source,thesmallerthespatialcoherencelength,andthereforethe
Proof of Claim 1. To evaluate the correlation, we note from highertheresolutionatwhichwecandiscriminatebetweenpaths
Equation(5)thattheFouriertransformofu(x,z,t)overthetime withdifferentpointsoforiginonthelightsource.
domain is the signal a eikz. Similarly the Fourier transform of
k ProofofClaim2.Aplanewavepropagatinginangleθis
u(x,z+τ,t)isa eikτeikz. UsingParseval’stheorem, theinner
k
productsinthetemporalandfrequencydomainsareequivalent, (cid:90)
u (x,z,t)= a ·exp(−ik(sin(θ)x+cos(θ)z−ct)) dk,
θ k,θ
(cid:28)(cid:90) (cid:29)
k
corr(x,τ)= |a |2·eikτdk . (12) (18)
k
whereweassumedη =1forsimplicityandwithoutlossofgener-
ality.Forsmallθ,sin(θ)≈θandcos(θ)≈1,therefore
Then,usingalsoEquation(6),wehave
(cid:90)
corr(x,τ)=(cid:90) e−(k2−∆k¯2k)2 ·eikτdk. (13) uθ(x,z,t)= ak,θ·exp(−ik(θx+z−ct)) dk. (19)
Thetimedelaybetweentheshiftedwaveu (x+ξ,z+τ,t)and
θ
Equation(13)is(uptoaconstant)theFouriertransformofaGaus- theoriginaloneu (x,z,t)istheprojectionoftheshift(ξ,τ)on
θ
sian with s.t.d ∆k. For a Gaussian centered at zero, the Fourier thepropagationdirection(sin(θ),cos(θ))≈(θ,1),
transformisaGaussianwithinverses.t.d. 1/∆ . AstheGaussian
k
isnotcenteredatzero,thisismultipliedbyasinusoidrepresenting u (x+ξ,z+τ,t)=u (x,z,t−(θξ+τ)/c). (20)
θ θ
2µm
AsmalladaptationofClaim1showsthat,foreveryθ, 0.9 3nm
0.8 10nm
corrθ(x,ξ,τ)=eik¯(θξ+τ)e−12(∆k(θξ+τ))2, (21) sity00..76 25nm
wheretheshiftbetweenthewavesisnowθξ+τ insteadofτ as en0.5
in Claim 1. If this shift is large, the interference contrast for ev- nt0.4
ery θ decays as a Gaussian with width equal to the temporal co- i0.3
herence length 1/∆ . However, even if the temporal coherence 0.2
k
lengthislargerelativetothespatialshift,sothatinEquation(21) 0.1
e−21(∆k(θξ+τ))2 ≈ 1,theobservableinterferencestilldecaysifξ −20 −10 0pathl1en0gth2d0iffer−en8c.0e−τ7(µ.5m−)7.0 −6.5−6.0
is larger than ∆ . To see this, we compute the total interference
c Figure4: Measuringtheopticalcoherencelength. Thegraphto
contrastbyintegratingEquation(21)overanglesinthesource,
theleftshowstheintensitymeasuredbythecamerainFigure2for
corr(x,ξ,τ)=(cid:90) eik¯(θξ+τ)e−12(∆k(θξ+τ))2dθ. (22) gdrifafeprhentot vthaelureisghotfsphaotwhlsenagcthlosdeif-fuepreonvceerτar=ang|deso−f2dµr|m, a.nWdhthene
θ thepathlengthdifferenceiszero,maximalconstructiveinterference
exists. Apathlengthshiftofλ/2resultsindestructiveinterference
ThoughwecancomputeEquation(22),toobtainasimplerexpres-
sionweapproximatetheGaussiantermbyitscentralvalue, (zerointensity). Asthepathlengthdifferenceincreases,correlation
betweenthe wavesis reducedand thefringe contrastdecays, be-
corr(x,ξ,τ)≈e−12(∆kτ)2(cid:90) eik¯(θξ+τ)dθ cTohmeidnigffearlemnotlsytczoerlooroedncpeloτtsexshceoewdsmtehaesutermempoenratslfcoorhderifefnerceenltesnpgetch-.
θ
tral bandwidths. Using a white light source with color filters of
=W∆c(ξ)eik¯(τ)G∆−1(τ). (23) spectralbandwidthof3,10and25nmproducestemporalcoher-
k
encelengthsof50,25,and10µm,respectively.Themeasurements
Thelastequalityfollowsfromthefactthatthetermeik¯θξisasinu- forthisfigureweretakenataresolutionof10nm.
soidvaryingwithθ,whoseintegralisasinc. Therefore,forlarge
ξ,integratingoverθcancompletelyeliminatethecorrelation.(cid:3)
bundleofmultiplepinholelasersources,orusingspeciallymanu-
3.3 CoherencePropertiesofDifferentSources factured,large-arealasersemiconductors(diodebarsandstacks).
There are many types of light sources with different coherence
Practical coherence values. To illustrate how these calculations
properties. Aswewillseeinthenextsection(e.g. Figure5),these
translate into practice, consider the coherence lengths of an LED
properties determine the type of transport decomposition that our
sourcecoupledwithacolorfilterandaperture.
interferometerprovides. Wediscussvarioussourceshere,limiting
ourattentiontothoseinthevisiblespectrum.
With respect to temporal coherence length, note that broadband
Lasersourcesareatoneextremeendofcoherence.Thesearehigh-
sources and color filters are often described by the width ∆ of
powersourcesthatcantypicallybemodelledasidealpoints(zero thevariationaroundthecentralwavelengthλ¯,ratherthanthevλari-
area), while also having very narrow, effectively monochromatic,
ation∆ aroundthewavenumber.UsingsimpleTaylorexpansion,
spectralwidths. Asaresult, theirtemporalandspatialcoherence k
wecanrelatethistothetemporalcoherencelengthas
lengths can reach a few meters, implying poor spatial and path-
length resolution. For this reason, lasers are not appropriate for λ¯2
micron-scalelighttransportdecompositions,butcanbeusedwhen L =1/∆ ≈ . (24)
micron-scaleresolutionisnotnecessary[Abramson1983]. c k 2π∆λ
Attheotherextremearelight-emittingdiode(LED)sources.These Imagine a broadband LED, emitting in the entire visible range
can be either “colored”, with a spectral bandwidth of a few tens [400,700] nm. Ignoring any ultraviolet and infrared energy, we
ofnanometersaroundtheircentralwavelength,orverybroadband, can set λ¯ = 550nm and ∆ = 150nm, resulting in coherence
λ
coveringtheentirevisiblespectrum. Necauseoftheirlowerpower lengthLc ≈ 0.3µm. Inpractice, thisveryfineresolutioncomes
density,LEDsourceshavelargeremittingareas,withwidthsrang- atthecostofveryreducedinterferencecontrast,duetochromatic
ingfromseveralmillimeterstoafewcentimeters.Asaresult,they aberration of optics. Alternatively, we can couple such a source
have micron-scale temporal and spatial coherence lengths. Their withacolorfilter,ordirectlyuseacoloredLEDsource. Theseal-
exactvaluescanbecontrolledbymodulatingthesourceoutputwith ternatives give spectral widths in the range of 1−50nm, which
anopticalcolorfilter,andplacinganaperturebetweenthesource correspond to temporal coherence lengths of a few microns. Ad-
and the collimating lens: A narrower spectral filter means longer ditionally, asmentionedinSection2, thesenarrowerbandsallow
temporalcoherence,andasmalleraperturemeanslargerspatialco- capturingmeasurementsatmultiplespectralchannels. InFigure4,
herencelength. Theincreaseincoherencecomesattheexpenseof we use the Michelson interferometer of Figure 2 to measure the
lowerlightoutput,andthereforelongerexposuretime. temporalcoherencelengthsofsuchacombination,usingthreedif-
ferentcolorfilters. Usingafilterofbandwidth25nm,wemeasure
Between the two extremes, there are two source categories. The atemporalcoherencelengthofabout10µm.Thisistheresolution
firstcategoryincludessuperluminescentdiode(SLD)sources,and weuseinmostofourexperimentsinthepaper.
supercontinuumlasersources.Likestandardlasers,thesearehigh-
powersourcesthathaveinfinitesimaleffectivearea,butatthesame
time,theyarepolychromaticlikestandardLEDs.SLDshaveband- Intermsofspatialcoherencelength,asindicatedbyEquation(17),
widths comparable to colored LEDs, and supercontinuum lasers thevalueswecanachieveinpracticeareafunctionoftheangular
spantheentirevisiblespectrum. Therefore,thesesourcescombine extent of the source. For visible wavelengths, assuming that we
themeter-sizespatialcoherencelengthsoflaserswiththemicron- cancollimatethesourceoutputtohaveanangularrangeof4◦,we
scaletemporalcoherencelengthsofLEDs. have ∆c ≈ 8µm. With an angular extent of 1◦ the coherence
lengthincreasesto∆ ≈ 32µm. Inanaturalenvironmentwhere
c
Finally,theoppositecombinationcorrespondstoamonochromatic illuminationarrivesfrommultipledirectionsoverthehemisphere,
areasource. Suchasourcecanbecreatedusingaspatiallydense coherencelengthcanapproachthewavelengthoflight.
4 Light Path Decomposition siblepaths,withatimedelayproportionaltothepathlength,
Equipped with analytical representations of temporal and spatial (cid:90) (cid:90) (cid:18) (cid:96)(o )(cid:19)
coherence,wearereadytointroduceamathematicalmodelofin- us,θ(x,t)= α(oγ)uin,θ x+ξ,t− cγ dξ
terferometricdecompositionoflighttransport.Thecoreofthissec- ξ {oγ}γ∈A(x+ξ,x)
tionisEquation(34),whichrelatesdiscreteinterferencemeasure- (cid:90) (cid:90) (cid:16) ς(cid:17)
ments to an underlying pathlength-resolved transport function—a = T(x+ξ,x,ς)uin,θ x+ξ,t− c dςdξ. (27)
continuous version of the pathlength-resolved transport matrix of ξ ς
Equation(3). Usingthismodel,were-interprettheoutputoftime-
The above equation relates the field to the function T(x ,x ,τ)
domain,full-fieldOCTasdecompositionbypathlengthτ(transient forthescene.Thetermu (cid:0)x+ξ,t− ς(cid:1)iscomplex-val1ued2and
imaging), andwedescribecombinationsofsourcesandreference in,θ c
includesthephasedifferencebetweentheintegratedcomponents,
mirrorsthatprovideothertypesofdecompositions.
duetothedifferentlengthsofpathsineachcomponent.
Throughoutthissection, weassumethatweimagewithacamera
focusedatdepthz ,asshowninFigure3. Aswerestrictthedis- Interference of reference and scattered fields. Having derived
o
cussiontopointsonthez = z plane,fromthispointforwardwe expressionsforthereferenceandscatteredfields,wenowcompute
o
simplifynotationbyomittingthezcoordinateoffields.Wedenote the interference contrast that can be measured by a camera. Ig-
byu (x,t)theincomingfieldarrivingatboththereferenceand noringfornowdiffractionblur,acamerafocusedatdepthzo will
targeitn,aθrmsattimetafterbeingemittedfromthesource. Thisis measurethesuperpositionofthetwofields,
asuperpositionofmonochromaticwaves,asinEquation(18). We (cid:90)
(dxen,zoote)bayndusti,mθ(ex,t.t)Stihmeifilaerlldy,swcaettdeerendotbeybtyheusr,cθe,nτe(xa,tts)ptahceerpeofeinr-t Ii,τ(x)= θ(cid:10)|ur,θ,τ(x,t)+us,θ(x,t)|2(cid:11)t dθ (28)
encefieldatpoint(x,zo)whenthereferencemirrorispositioned =Is(x)+Ir,τ(x)
atdepthz +τ/2.Webeginbyderivingexpressionsthatrelatethe
o (cid:26)(cid:90) (cid:27)
referenceandscatteredfieldstotheinputfield. +2Re (cid:104)u (x,t)∗·u (x,t)(cid:105) dθ , (29)
r,θ,τ s,θ t
Referencefield. Thereferencefieldu relatestou through θ
r,θ,τ in,θ
aspatialtransformationf(x)andatemporalshift, wherebyI (x),I (x)wedenotetheintensitiesofthescattered
s r,τ
andreferencefields,respectively,
u (x,t)=u (f(x),t−τ/c). (25)
r,θ,τ in,θ
(cid:90) (cid:90)
InthestandardMichelsoninterferometerofFigure2,wehavesim- Is(x)= (cid:10)|us,θ(x,t)|2(cid:11)t dθ, Ir,τ(x)= (cid:10)|ur,θ,τ(x,t)|2(cid:11)t dθ.
plyf(x)=x.Wecanproducemoregeneraltransformationsf(x) θ θ
(30)
withothermirrorconfigurations(Figure6). Theformoff(x)will
TheintensitycomponentsI andI correspondtotheimagesthe
determinethenon-zeroentriesofthematrixMfromEquation(3). s r,τ
camerawouldmeasureifweblockedthemirrororscenearmsofthe
Scattered field. The scattered field u relates to the input field setup, respectively, and measured the other. By subtracting these
s,θ
throughamoregeneraltransformation,whichisafunctionofthe components,weremainwiththeinterferencecontrast,equaltothe
scene’s geometric, refractive, and material properties. To derive realpartofthecorrelationofthereferenceandscatteredfields,
its form, we introduce some additional notation. We denote by
(cid:90)
{(xoγ,}zγ∈)At(ox1p,xo2in)tth(xes,ezto)fa(slleepoFsisgiubrleep3afthosritnhethceassecexnef=romxp),oibnyt Cτ(x)= (cid:104)ur,θ,τ(x,t)∗·us,θ(x,t)(cid:105)t dθ. (31)
1 o 2 o 1 2 θ
(cid:96)(o )theopticallengthofeachpath,andbyα(o )theenergyloss
γ γ
alongthepath—thislossisafunctionofthepropagationdistance Cτ(x)isacorrelationsignal, analogoustocorr(x,ξ,τ)inEqua-
andvolumetricabsorptionalongthepath.Dependingonthescene, tion(15). Wedenoteitbyadistinctsymbolbecauseofitsimpor-
suchpathscanbeofmanyforms: directpathstravelingfromthe tanceinourapplication. Usingthecoherencepropertiesforcorre-
sourcetotheobjectandthenthecamera,multiplereflections,mul- lationsignalsfromSection3,wecanprovethefollowingclaim.
tiple scattering, and so on. In the simple case that the scene is a
Claim3. Thecorrelationofthereferenceandscatteredfieldsis
perfectplanarmirror,therearepathswithnon-zeroabsorptiononly
if the start and end point are the same, x1 = x2, and for every (cid:90)
x(x11o,nzloy).aFsoinrgalepastuhchinpaaitrh,wexhiesrtes—thteheredfriarecctitvpeaitnhdferxomη ≈(x11,,z(cid:96)o()oγto) Cτ(x)∝ ξW∆c(ξ)
isequaltothegeometricpathlength;inmediawithlargerrefractive (cid:90)
indices,thepathlengthalsofoldsintheopticalpathdelay.Finally, T(f(x)+ξ,x,ς)eik¯(ς−τ)G∆−1(ς−τ)dςdξ. (32)
we denote by T(x ,x ,τ) the complex scalar resulting from the ς k
1 2
integrationofallpathsfromx tox withopticallengthτ,
1 2
(cid:90) Proof. Wedenoteξ=x+ζ−f(x),implyingx+ζ =f(x)+ξ.
T(x1,x2,τ)= α(oγ). (26) SubstitutingEquations(25)and(27)inEquation(31),weobtain
{oγ|(cid:96)(oγ)=τ}γ∈A(x1,x2) (cid:90) (cid:90)
C (x)= T(f(x)+ξ,x,ς)
WecallT(x ,x ,τ)thepathlength-resolvedlighttransportfunc- τ
1 2
tion. Whenthesceneisaplanarmirror, T(x ,x ,τ) = δ(x − ξ ζ
x ,τ −d ), where d is the depth of the sce1ne2mirror. By d1ef- (cid:90)(cid:68) (cid:16) τ(cid:17)∗ (cid:16) ς(cid:17)(cid:69)
2 s s u f(x),t− u f(x)+ξ,t− dθdζdξ (33)
inition,weobservethatthefunctionT(x1,x2,τ)isthecomplex- θ in,θ c in,θ c t
valuedcontinuousversionofthepathlength-resolvedlighttransport
matrix,Tτ,introducedinSection1. Inpartsofthefollowingdis- UsingClaim2,wegetEquation(32).(cid:3)
cussion,wewillbeemployingthelighttransportmatrixdiscretiza-
tionforvisualizationandintuition,withtheunderstandingthatwe TheaboveclaimshowsthatthecorrelationCτ(x)isequaltothe
continuetousethecontinuousversionT(x ,x ,τ). pathlength-resolvedlighttransportfunctionT(f(x),x,τ)fromen-
1 2
trancepoint(f(x),z )toexitpoint(x,z ),blurredinspacewith
o o
We can now use the above functions to derive the scattered field a sinc of width proportional to the spatial coherence length, and
u (x,z ).Thisisthesuperpositionofcontributionsalongallpos- blurredinthepathlengthdimensionwithaGaussianofwidthpro-
s,θ o
SLD mirror diffuser SLD or a supercontinuum laser. Capturing an entire transient se-
spatially (cid:51) quencecorrespondstosweepingoverτ values,whichcanbedone
o
coherent (τ) bydenselyscanningthereferencearmtodifferentpositions. This
temporally (cid:55) hardwarecombinationandcaptureprocedureisexactlyequivalent
coherent toconventionalfull-field,time-domainOCT.
(a)pathlengthdecomposition
Spatial probing. Referring again to Equation (3), another form
laserbundle
of decomposition corresponds to setting, w(τ) = 1 everywhere,
spatially (cid:55)
andusingaspatialprobingpatternMwithnon-zeroentriesonly
coherent (cid:80)τ (τ) for the paths we want to preserve. This requires a source with a
temporally (cid:51) near-infinitetemporalcoherencelengthandaverysmallspatialco-
coherent herencelength. FromSection3.3,suchasourcecanbeproduced
(b)spatialprobing using a laser bundle. Additionally, we require a reference mirror
LED configurationwhoseshapeinducesthespatialtransformationf(x)
spatially (cid:55) correspondingtothedesiredpatternM.Usingaplanarmirrorasin
coherent (τ) thestandardMichelsoninterferometercorrespondstoamatrixM
temporally (cid:55) thatisnon-zeroonlyinitsmaindiagonal. Wediscussalternatives
coherent laterinthesection. Duetotheinfinitetemporalcoherencelength,
thiscasedoesnotrequireascanofthereferencearm. Wemention
(c)pathlengthdecomposition+spatialprobing
spatialprobinghereforcompleteness,butnotethat,forthistypeof
Figure5: ByequippingtheMichelsoninterferometerwithsources decomposition,interferometrydoesnotofferanyadvantagescom-
havingdifferentcoherenceproperties,wecancapturedifferentlight pared to the other, easier to implement and use, techniques dis-
path decompositions. (a) Pathlength decomposition separates all cussed in Section 2. For this reason, we did not implement this
pathsinascene(M = 1),regardlessoftheirendpointlocations, case,andonlyusedspatialprobingcombinedwithpathlengthde-
intermsoftheiropticallengthτ. Itseparatesdirectpaths(blue), composition,asdiscussedinthenextparagraph.
reflections (purple), retroreflections (orange), and all scattering
Pathlengthdecompositionandspatialprobing.Thelastcasede-
(greens). (b) Spatial probing separates paths with certain start
composeslighttransportsimultaneouslyintermsofpathlengthand
pointlocationswhileignoringtheirlengths(byeffectivelymeasur-
spatialcorrespondences. Assuggestedbytheprevioustwocases,
ingsummationsoverpathlength). ForadiagonalM,itseparates we can do this using a source with temporal and spatial coher-
direct, retroreflection, and backscattering (dark green) paths. (c) encelengthsthatarebothsmall.AnLEDsourceistheappropriate
Finally, itisalsopossibletocombinetheprevioustwocasesand choiceasdescribedinSection3.3.Asisneededforpurepathlength
separatepathsintermsofbothpathlengthandendpointlocations. decomposition,capturerequiresadensescanofthereferencearm
tomanydifferentpositions.
portionaltothetemporalcoherencelengthofthesource.Byrecord-
ingtheintensityI (x)ofEquation(28)andsubtractingtheindi- Creatingprobingpatterns. Inthelasttwodecompositiontypes,
i,τ
vidualintensitiesofthesourceandreferencefields,wehave replacingtheplanarmirrorMrinthereferencearmoftheMichel-
son interferometer with different mirror shapes induces different
Tm(f(x),x,τ)(cid:44)|I (x)−(I (x)+I (x))|2 (34) correspondences f(x), that can be used to probe off-diagonal el-
i,τ s r,τ
ementsofthelighttransportmatrix. Forinstance, theright-angle
=|2Re(Cτ(x))|2 pairofmirrorsinFigure6(b)inducesaverticalflip,orf(x)=−x.
=(cid:12)(cid:12)(cid:12)2Re(cid:16)T(f(x),x,τ)∗W∆c ∗(eik¯τG∆−1(τ))(cid:17)(cid:12)(cid:12)(cid:12)2. Tpohritsmcoartrreisxp,owndhsicthoimncelausduerisngtwthoe-baonutni-cdeiapgaotnhaslaosftshheowlignhitntrFanigs--
k ure 6(b). Shifting such a pair of mirrors vertically in the figure
Thisequationisourmainresult. Itshowsthatusinginterference, measuresarbitraryanti-diagonals,asinFigure6(c);andcombina-
we can measure blurred samples of the pathlength-resolved light tionsofparallelmirrorsandorthogonalmirrorpairssamplediffer-
transportfunction. Thespatialresolutioniscontrolledbythesize entblocksofthelighttransportmatrixindifferentways,asshown
of the spatial blur kernel, which is in turn controlled by the illu- inFigure6(d). NotethattheprobingpatternMproducedbythese
minationangleofthelightsource. Similarly,thepathlengthreso- mirror configurations will be different if we change our setup so
lutioniscontrolledbythesizeoftheblurkernelinthepathlength thatthecameraandsourcearenolongerco-axial.
dimension, which in turn is controlled by the spectral bandwidth
of the light source. The exact correspondence (f(x),x) between 5 Pipeline and Design Considerations
entranceandexitpointswherewesamplethelighttransportfunc-
tion depends on the transformation implemented at the reference Havingpresentedoursetupandmethodology,wenowdiscussthe
armofthesetup.Finally,thepathlengthwherewesamplethelight computational post-processing and other practical considerations
transportfunctiondependsonthedepthofthereferencearm. foroptimizingtheperformanceofoursetupforgraphicsapplica-
tions. Throughoutthissection,weconsideratoysceneconsisting
4.1 DecompositionTypes ofatilteddiffuseplane,asshowninFigure7(a).
BasedonourmodelinEquation(34),theinterferometerofFigure2 Capture and computational post-processing. Measuring the
canbeequippedwithvarioussourcesandmirrorconfigurationsto scene and mirror images I and I separately is impractical, as
s r,τ
producedifferenttypesoflighttransportdecompositions. Thefol- it doubles acquisition time. In practice, since the reference arm
lowingdiscussionissummarizedinFigures5and6. contains a mirror configuration, the intensity I does not vary
r,τ
muchwithτ andcanbeestimatedbysmoothingtheτ dimension.
Pathlengthdecomposition.Thefirstcaseisequivalenttotransient
Forthisreason,weperformasinglescanofthescene,capturinga
imaging,andcorrespondstosettingM=1andw(τ)=δ(τ−τ )
o setofframescorrespondingtoI (x)forvariousvaluesofτ, as
inEquation(3).Fortheformer,werequireasourcewithverylarge, i,τ
shown in Figure 7(a)-(b). Scanning a volume of 1cm using this
essentiallyinfinite,spatialcoherencelength,meaningthatthespa-
processtypicallyrequirescapturing10000images(oneimageper
tialblurkernelW inEquation(34)isequalto1everywhere.For
∆c micron-translationofthereferencearm),aprocesswhichcantake
thelatter, weneedasourcewithaveryshorttemporalcoherence
uptoseveralhoursforscenesrequiringlongexposuretimes.Aswe
length,makingthepathlengthblurkernelG inEquation(34)
∆−1 mentioninSection6,forscenesthatrequiremultiplespectralchan-
k
verynarrow. AsdiscussedinSection3.3,thisisachievedusingan nelsorhigh-dynamicrange,weneedtorepeatthescanningprocess
mirror diffuser mirror diffuser
tilteddiffuse
plane
τ
N
τ1 ... τN ...
capturedframes τ1
(a)diagonalprobing (b)anti-diagonalprobing
τ ... τ
1 N
(c)sub-anti-diagonalprobing (d)otherprobingpatterns processedframes
Figure6:ReferencemirrorconfigurationsinaMichelsoninterfer- Figure7: Captureandcomputationalpipeline. Weconsideratoy
ometerinducedifferentspatialprobingpatterns.Assumingco-axial sceneconsistingofadiffuseplane,tiltedrelativetotheopticalaxis
sourceandcamera: (a)Aplanarmirrorprobesthemaindiagonal of the camera and illumination, as shown in the upper left. By
of the light transport matrix, corresponding to direct, retroreflec- translating the reference mirror to different positions, we capture
tion,andbackscatteringpaths.(b)Aright-anglemirrorpairprobes a sequence of frames corresponding to the intensity Ii,τ(x) mea-
theanti-diagonal, whichincludes2n-bouncereflectionpaths. (c) suredfordifferentvaluesτ. ByprocessingtheframesusingEqua-
Offsetting a right-angle pair allows probing away from the anti- tion(36),weseparatetheinterferencecomponentTˆm(f(x),x,τ),
diagonal. (d)Differentcombinationsofthesebasicconfigurations resultingintheimagesatthebottom. Asthetargetplaneistilted,
probeblocksofthelighttransportmatrixindifferentways. different points have different depths and therefore interfere with
the target arm at different pathlength values τ, corresponding to
severaltimes,resultinginhundredsofthousandsofimages. thescanlinesshownintheprocessedframes.
Followingthecapturesection,inpost-processingweapproximate
theinterference-freecomponentas theenergyofthescatteredfieldattheopticallengthτ ofinterest
Tm(f(x),x,τ). Becauseindiffusescenesonlyasmallportionof
Iˆ(x)+Iˆ(x)≈median (I (x)). (35) the energy in the target arm returns to the camera, to achieve the
s r τ i,τ
above ratio, we need to attenuate the very high amplitude of the
Furthermore, in real scenes the interference signal Tm of Equa- referencearm. Inmostscenes,theintensityatasinglepointisthe
tion34involvesahigh-variationpseudo-randomspecklenoise. To resultofcontributionsfrommultiplepaths,makingitimpossibleto
eliminatethesespeckleartifacts,weblurourmeasurementsovera a-priorimatchtheenergyofthepathlengthτ componentonly. In-
smallspatialwindowwithafilterg(x). Ourfinalestimateofthe stead,weattenuatethereferencearmsothatitsintensitymatches
opticalpathlengthdecompositionfunctionis,then,computedas thetotalintensity.AsdiscussedinSection6,wedothisusingeither
cross-polarization,orneutraldensityfilters.
Tˆm(f(x),x,τ)=Tm(f(x),x,τ)∗g(x)
Opticsblur. TheidealanalysisoftheSection4doesnottakeinto
=|I (x)−Iˆ(x)−Iˆ(x)|2∗g(x). (36)
i,τ s r accountblurintroducedbytheopticsoftheinterferometrysetup,
includingdiffractionandsensorblur. Theseblurprocessescanre-
Thisblurringisoneofthefactorsthatresultinareductioninthe
duceinterferencecontrast,unlessthevariousapertureandmagni-
spatialresolutionofourmeasurements.Wediscussspecklenoisein
ficationparametersarechosenappropriately. Wediscusstheseef-
moredetailinthesupplementarymaterial. Intheremainingofthis
fectsbelow,andingreaterdetailinthesupplementarymaterial.
section,wediscussotherfactors,aswellaspracticestomaximize
thesignal-to-noiseratioofthemeasurements(36).
Weconsiderasbeforeacamerafocusedatz ,asshowninFigure3.
o
Contrast. Assumewewanttocapturethepathlengthcomponent Fromdiffractiontheory[Goodman1968],wehavethatinsteadof
whichcorrespondstoopticallengthτ andthatthereferencemirror thefieldu (x,t)arrivingatpoint(x,z ),thecamerameasuresthe
θ o
ispositionedatthecorrespondingdepth.Wewilldefinethecontrast fieldblurredbythecameradiffractionblurkernel,W ∗u (x,t),
asthespatiallyaveragedmagnitudeoftheinterference(correlation) wherethediffractionblurwidth∆ =λ¯/(2Φ)isinv∆erΦselypθropor-
Φ
component,overtheintensitycomponents, tionaltotheacceptanceangleΦofthecamera.
(cid:112)
Contrast= Ex[Tm(f(x),x,τ)]. (37) Theintensitymeasuredbythecameraistheintensityoftheblurred
E [I (x)+I (x)] field,averagedovertheareaofapixelonthesensor,
x s r,τ
Toavoidsaturation,wesettheexposuretimesothattheaveraged (cid:90) Θ/2
intensities without interference fall in the middle of the sensor’s I(x)=Π∆x ∗ (cid:10)|W∆Φ ∗uθ(x,t)|2(cid:11)t dθ, (38)
dynamicrange,E [I (x)+I (x)]≈0.5.Then,Equation(37)is −Θ/2
x s r,τ
simplytheaveragedmagnitudeoftheinterferencecontrast.
whereΠ isthepixel-sizedrectangularfunction.
∆x
AshortcalculationshowsthatthecontrastinEquation(37)ismax-
imized when the energy of the reference field I (x) is equal to Asaresultoftheabovetwoprocesses,weproveinthesupplemen-
r,τ
0.45 optimized 0.45 optimized
2×A
(b) (c) (d)(e) (g) (h) 0.40 4×As 0.40 2×As
0.35 2×ps 0.35 4×As 2×Ac
4×p
As (f) ntrast00..2350 24××AAcc ntrast00..2350 42××Apc
(j) o0.20 o0.20
(a) Θ/2 (i)c c
Ac 0.15 0.15 4×p
(k) 0.10 0.10
0.05 0.05
Φ/2 p
(l) −37.85−37.84−37.83−37.82−37.81 0 0.2 0.4 0.6 0.8 1
pathlengthτ (µm) exposuretime
Figure8: OptimizingtheMichelsoninterferometerofFigure2. Theschematictotheleftdescribesourimplementation. Weconsiderthe
effectofthefollowingimagingparameters: sourceapertureareaA ,cameraapertureareaA ,andpixelpitchp. Theoptimizedsettings
s c
correspondtocameraF-numberoff/32,sourceaperturesothatthesource’sangularextentisapproximatelyhalfthatofthecamera,and
pixelpitchof4µm.Differentconfigurationsareproducedbymodifyingonlyoneparameter.Foreachconfiguration,themiddlegraphshows
themeasuredcontrastforthetilteddiffuseplanesceneofFigure7;andtherightgraphcomparescontrastandtotalexposuretime.
tarymaterialthatmeasurementsfromEquation(34)taketheform, toautomaticallychangebetweenfilterscenteredatdifferentwave-
lengths:625nm(red),525nm(green),and450nm(blue).
Tm(f(x),x,τ)=|2Re(β∗C (x))|2
τ Additionally,weusethreepolarizers: atthesource(d),camera(j),
=(cid:12)(cid:12)2Re(cid:16)β∗W ∗T(f(x),x,τ)∗(eik¯τG (τ))(cid:17)(cid:12)(cid:12)2, andreferencearm(g).Thepolarizersonthesourceandcameraare
(cid:12) ∆c ∆−1 (cid:12) rotatedtobeeitherparallelorcrossedwitheachother,depending
k
(39) onwhichpolarizationcomponentofthescene’stransportistobe
measured.Thepolarizeronthereferencearm(g)isrotatedrelative
whereβ =Π∆x ∗W∆Φ. totheothertwoinordertoattenuatethereferencebeambydifferent
amounts. ThisenablesHDRimaging: wescanexposurebrackets
Itisimportanttoobservethattheadditionalblurβactsonthecom-
where,foreachsettingofthecamera’sexposuretime,thereference
plex valued correlation component, rather than on powers (inten-
polarizer (g) is adjusted so that the reference beam’s intensity is
sities). Therefore, ifthecorrelationsignalC (x)containsspatial
τ at roughly one-quarter of the sensor’s dynamic range (for a total
featureswhosefrequencyishigherthanthewidths∆ ,∆ ofthe
x Φ intensityatroughlyhalfthedynamicrangewhenimagedtogether
blurkernels,theaveragedpowerofmeasurementsignalisreduced.
with properly exposed parts of the target scene). When we need
FromEquation32,weknowthatthesignalC (x)haslimitedspa-
τ tomeasurebothpolarizationcomponents,weremovepolarizers(d)
tialresolution∆ ,asaresultofthelimitedspatialcoherenceofthe
c and(j),andreplace(g)withasequenceofneutraldensityfilters.
lightsource.Therefore,toachievegoodcontrast,weneedtoselect
theimagingparameterssuchthattheopticsblurremainsbelowthe Wediscussourimplementationindetailinthesupplementaryma-
signalresolution,thatis,∆x <∆cand∆Φ <∆c. terial. Intherestofthissection, weusethisassemblytoanalyze
transport in various scenes and for various applications. We first
Wenotethatrequesting∆ < ∆ isequivalenttorequestingthat
Φ c discussexperimentsusinganLEDsource(Figure5(c)).
theangularextentofthesourceissmallerthantheacceptanceangle
ofthecamera,Θ<Φ.Ingeneral,usingawidesourceisdesirable, Depthscanning.FollowingSection4.1,weuseanLEDsourceand
asitincreasestheoverallpoweroftheilluminationandtherefore aplanarreferencemirrortomeasurethepathlengthdecomposition
reducesexposuretime.However,increasingthesourceangularex- ofthediagonalcomponentTm(x,x,τ).Asourcameraandsource
tentwithoutadjustingtheimagingparametersaccordinglywillalso areco-axial,wethenobtainthedepthD(x)atpixelxas
reducethecontrastofthecorrelationsignal.Torelatethesedirectly
toimagingparameters,wedenotebyΨthenumericalapertureof D(x)=min{τ :Tm(x,x,τ)>0}. (40)
thecameralens.Ashortcalculationshowsthatacameraimagingat
magnificationmacceptslightuptoamaximalangleΦ= Ψ . Figure 9 shows scans of a few different objects. Each scan mea-
1−1/m sures a volume of 15mm at a step size of 1µm, and requires
InFigure8,wequantifytheaboveobservationsbymeasuringthe roughlyeighthoursofcapturetime. Theeffectivedepthresolution
detectioncontrastforthetilteddiffuseplanesceneofFigure7un- is10µm. Wecanmeasuredepthinanumberofchallengingsitu-
derdifferentimagingconfigurations.Wevarythesizeofthesource ations,suchaswhenmaterialsexhibitsubstantialscattering(soap),
and camera aperture, and the size of the camera pixel—which is semi-transparency(gummybear),andstrongcausticsandspecular
analogoustochangingthecameramagnification. Weobservethat ordiffuseinterreflections(cupandpasta).Theveryfinepathlength
settingimagingparameterssub-optimallycanresultinlossofmore resolutionofoursetuptranslatesdirectlyintosharpdepthresolu-
thanhalfthedetectioncontrast. Ontheotherhand,optimizingde- tion,whichallowscapturingdetailssuchasthehairtextureonthe
tectioncontrastcomesatthecostofincreasedexposuretime. coinandthefinetextureonthegnocchisurface.
Figure 10 shows depth maps obtained for two more scenes. The
6 Implementation and Experiments
stagestepsizeandeffectiveresolutionarethesameasbefore,but
the volume and capture time increase to roughly 25mm and 12
Figure8showsaschematicofouropticalsetup.Weselectsources
hours,respectively,tocapturethelongerreflectionpaths. Thefirst
andmirrorscorrespondingtoconfigurations(a)and(c)ofFigure5,
sceneincludessharpspecularinterreflections, with asinglechess
and(a)and(b)ofFigure6, respectively. Theexactconfiguration
piecebetweenanangleddiffusewall(right)andamirror(left). In
dependsontheapplication,asdescribedbelow.
the acquired depth map, the mirror reflections are interpreted as
OurschematichastwoadditionscomparedtoFigure2. First, we realobjectsatagreaterdepth.Despitetheindirectreflectionbythe
modulatetheoutputofthesourcesourceusingcolorfiltersofband- mirror, the depth map of the chess piece is accurately recovered.
width25nm(e),producingatemporalcoherencelengthofabout The second scene shows a translucent gummy bear between two
10µm(Figure4). ToenableRGBimaging,weuseafilterwheel diffusewalls.ComparedtotherightmostcolumnofFigure9,where
Figure9:Scanningdepth.Foreachscenefromtoptobottom:anunprocessedHDRframe;scanneddepth,visualizedusingpseudo-colorand
aperiodicchangeinluminanceforacontouringeffect;andrenderingofscanned3Dmesh.Allscansareatadepth/pathlengthresolutionof
10µm.Scenesfromlefttoright:USquarter,gnocchi,soapcarving,plastictoycup,drypasta,andgummybear.
HDRframe directcomponent globalcomponent
Figure11:Direct-globalseparationforaclose-upofastrawberry.
theexactdirectcomponent,anditisnecessarytoalsodopathlength
decomposition. WithreferencetoFigure1,spatialprobing(diag-
onalofthelighttransportmatrix)capturesthevariousbluepaths,
whichincludedirectpathsandretroreflections, butalsobackscat-
Figure10: Depthscansofcomplexscenes. Fromlefttoright: an
tering. Usingpathlengthdecompositiontoselectonlytheshortest
unprocessedHDRframe;scanneddepth,visualizedasinFigure9;
pathadditionallyremovesthebackscatteringcontributions.
andrenderingof3Dmesh.Toprow:Chessknightbetweenadiffuse
walltotherightandamirrortotheleft. Themirror-reflectionis Figure 11 shows such a decomposition for a close-up of a straw-
interpretedasarealobjectatagreaterdepth.Bottomrow:Gummy berry.Weproducedthisseparationfrom12scansofthestrawberry
bearbetweentwodiffusewalls. Thescenegeometryisrecovered (4exposures×3spectralchannels).Eachscanwasperformedata
accurately, despite the conflation of paths from the diffuse walls stepsizeof10µm,meaningthatweslightlyundersampledinthe
behindthenear-transparentgummybear. pathlengthdomainrelativetothesource’spathlengthresolutionof
10µm.Additionally,weused8×8pixelbinningonthecamerato
the gummy bear is surrounded by air, here the combination of a speed-upframestreaming, andlargesourceandcameraapertures
diffusewallbehindanear-transparentgummybearinducesmany todecreaseexposuretime. Eventhoughtheseimagingsettingsare
distinctpathsofsimilarintensityateachcamerapixel.Despitethis suboptimal(seediscussionofFigure8),theywerenecessarytore-
conflation,thescenegeometryisstillrecoveredaccurately. duce the total scanning time to a little below two hours, thereby
accommodatingthelimitedshelflifeofthestrawberry.
Separation of direct and global components. We can use the
same source and mirror combination as before to also do direct- Weobservethatthedirectcomponentisessentiallyachromaticand
global separation. Since the source and camera are coaxial, the captures features characteristic of direct illumination, such as the
direct component at each pixel is the energy of the shortest path specular highlights and high-frequency surface structure. The in-
thatcontributestothediagonalofthelighttransportmatrix.Thus, directcomponentinthisimageisprimarilyvolumetricscattering,
andweseethatitaccountsfortheredcolorofthestrawberryanda
I (x)=Tm(x,x,D(x)), (41) largepartofthecolorationofthestrawberryseeds.
direct
whereD(x)isthedepthfromEquation(40). Theglobalcompo- Scattering. AsdiscussedinSection4.1,wecanalsoobtainade-
nentsumsthecontributionsfromallotherpathsand,givenanimage compositionofscatteringpathsthatbeginandendatthesamelo-
I(x)capturedusingaconventionalcamera,canbeobtainedas cationsonthesurfaceofascatteringvolume. Todemonstratethis,
Figure12showstheresultsfromscanningascenewithmaterialsof
I (x)=I(x)−Tm(x,x,D(x)). (42) differentscatteringproperties,suchasametal,adiffusedwallanda
global
highly-scatteringzirconiacoating. Figure12(c)showsaprocessed
Wenotethatspatialprobingaloneisnotsufficientforcomputing framefromthescan,wheretheintensitytailsatthezirconiacoat-
Description:build a prototype to measure and visualize three-dimensional shape, direct and cal interferometry to produce high-fidelity light transport decom- positions. Spatial decomposition methods induce different matrices M while keeping w(τ) trance point (f(x), zo) to exit point (x, zo), blurred in spa