Deconvolution using Fourier transform phase, ℓ1 and ℓ2 balls, and filtered variation

Contents lists available atScienceDirect

Digital

Signal

Processing

www.elsevier.com/locate/dsp

Deconvolution

using

Fourier

transform

phase,



₁

and



₂

balls,

and

filtered

variation

Onur Yorulmaz

a

,

∗

,

A.

Enis Cetin

b

,

1

a_{Dept.ofElectricalandElectronicsEngineering,BilkentUniversity,06800Bilkent,Ankara,Turkey} b_{Dept.ofElectricalandComputerEngineering,UniversityofIllinoisatChicago,Chicago,IL,USA}

a

r

t

i

c

l

e

i

n

f

o

a

b

s

t

r

a

c

t

Articlehistory:

Availableonline3January2018

Keywords:

Deconvolution Fouriertransformphase

1ball

Projectionontoconvexsets

Inthisarticle,wepresentadeconvolutionsoftwarebasedonconvexsetsconstructedfromthephaseof theFourierTransform,bounded2energyand1energyofagivenimage.Theiterativedeconvolution algorithmisbasedonthemethodofprojectionsontoconvexsets.Anotherfeatureofthemethodisthat it canincorporateanapproximatetotalvariationboundcalledfilteredvariationboundontheiterative deconvolutionalgorithm.Themainpurposeofthisarticleistointroducetheopensourcesoftwarecalled projDeconvv2.

©2018ElsevierInc.Allrightsreserved.

1. Motivationandsignificance

Deconvolutionistheprocessofinvertingtheeffectsoffiltering that reduces the image quality – mostly by blurring – in many image processing applications. The source of the blurriness may vary form camera motion to optical characteristics of the image capturingequipment. Thereforeit may be crucialin manyimage processingproblemstodeconvolve an inputimage beforefurther processing.

In this article we present a software that uses Projections ontoConvexSets (POCS)basedmethodinordertocorrecthighly blurred out of focus images. It may not be possible to estimate the parameters of blurring process in out of focus microscopic andMagneticParticleImaging(MPI)images.Insuchhighlyout-of focus images,well-known deconvolutionalgorithms are not very efficient. In these imaging problems it may not be possible to estimate a point spread function (psf) which is crucial for most deconvolutionalgorithms.However,insomecasesitmaybe possi-bletoassumethatthepsfissymmetricwithrespecttotheorigin. Asa resultitis possibleto estimatethe phaseofthe input from theobservedimageifthepsf issymmetric.Oursoftwareexploits thesymmetry characteristicsofthepsf andestimatesthe Fourier Transformphasesfromtheblurryimage. Softwarethen uses iter-ativePOCSmethodonFourierphase,



1 and



2 balls,andFiltered

Variationinordertoperformdeconvolution.

*

Correspondingauthor.

E-mailaddress:yorulmazonur@gmail.com(O. Yorulmaz).

1 _{A.EnisCetinisonleavefromBilkentUniversity.}

POCS baseddeconvolution was first developed by Trusseland Civanlar [19]. The method relies on iterative projections onto known convex properties of the image in spacial/frequency do-mains. Let the observed image y be the blurred version of the original image x0 andh be the blurringfunction. In many

prob-lems y

[

n1

,

n2

]

isalsocorrupted bynoise.Foragivenimage pixel

[

n1

,

n2

]

wedefineahyperplaneasfollow:

Cn1,n2

=

⎧

⎨

⎩

x|y[n1

,

n2

] =



k1,k2 h[k1

,

k2

]x[n

1

−

k1

,

n2

−

k2

]

⎫

⎬

⎭

(1)

Thereforethesolutionofx0

[

n1

,

n2

]

mustbeintheintersection

ofthesehyperplanes:

x0

∈ ∩

n1,n2Cn1,n2 (2)

Projection onto Cn1,n2 is essentially equivalent to making

or-thogonalprojectionsonto hyperplanesbecausetheconvolution is alinearoperation.Letxibethecurrentiterateoftheiterative de-convolutionalgorithm.Itsorthogonalprojectionxi+1 ontoCn1,n2 is

givenby:

xi+1

=

xi

+ λ

y[n1

,

n2

] − (h

∗

x

)

[n

1

,

n2

]

||

h

||

2 h

,

i

=

1,2,3, ... (3)

where

λ

=

1 for orthogonal projection and x1 is an estimate of

x0.POCStheoryallows 0

< λ

<

2 forconvergence[22,15].Eq.(3)

abusesthenotationalittlebitbecausesizeoftheimage x andthe blurh maybedifferent.Theh vectorshouldbepaddedwithzeros beforeadditionoperation.

https://doi.org/10.1016/j.dsp.2017.11.004

But making successive orthogonal projections onto the sets Cn1,n2 maynot be sufficient to reconstruct theoriginal image x0

becauseof the ill-posed structure of the problem andthe noise. Thereforeotherclosedandconvexsetsrestrictingthesolutioninto afeasiblesetmaybenecessary.TrusselandCivanlarused



2norm

based regularizing sets in their POCS based deconvolution algo-rithm.Inouralgorithm,weusethephaseoftheFourierTransform and



1-ball,filtered variationbased setsinaddition to



2-ball.In

whatfollowswedescribeotherclosedandconvexsetsthatcanbe usedindeconvolutionproblems.

1.1.Fouriertransformphase

Inmanydeconvolutionproblems[19,2]theblurringfunctionis symmetricwithrespecttoorigin,i.e.,

h[n1

,

n2

] =

h[−n1

,

−n

2

]

(4)

For example, this condition is satisfied by Gaussian blurs, disk shaped blursand some motion blurkernels [5]. Such kernels do notdistortthe Fourierphase oftheinputimage.Thismeans that thephaseoftheobservedimageandtheoriginalimagearerelated andthephaseoftheoriginalimagecanbeestimatedfromthe ob-servedimage. Any iterativedeconvolutionalgorithmcan take ad-vantageofthisrelationbyperformingorthogonalprojectionsonto thephasesetduringiterations[18,21].Intheabsenceofnoise:

Y

(

ejω

)

=

H

(

ejω

)

X0

(

ejω

)

(5)

Ifthe bluris symmetric H

(

ejω

₎

_{is real.When H}

₍

_ejω

₎

_>

_{0 for} some

ω

1

,

ω

2,

 Y

(

ejω

)

=

 X0

(

ejω

).

(6)

WhenH

(

ejω

₎

_<

_{0 forsome}

_ω

1

,

ω

2,

 Y

(

ejω

)

=

 X0

(

ejω

)

+

π

.

(7)

Asaresultwecandeterminethephaseof X0 fromY

(

ejω

)

.The

setofimageswithagivenFourierTransformphaseisaclosedand convexset [22,12]. Therefore we introduce the following set for deconvolutionproblems Cφ

=

x

|

 X

(

ejω

)

=

 X0

(

ejω

)

(8) whichisthesetofimageswhoseFouriertransformphaseisequal toagivenphase X0

(

ejω

)

.

ProjectionontoCφisobtainedintheFourierdomain.Letxφibe theprojectionofxi.TheFouriertransform Xφi ofxφi isgivenby Xφi

(

ejω

)

= |

Xi

(

ejω

)

|

ej X0(ω) (9) wherethephaseofXiisreplacedbythegivenphase X0

(

ω

)

.The

spatialdomainimagexφi

[

n1

,

n2

]

willbe F−1

[

Xφi

(

ejω

)

]

where F−1 istheinverseFourierTransformoperation.Recently,thesubspace orsparsity constraintswere imposed to reduce the search space

[9,1].The phase set is alsoa subspace whichreduces the search space.

1.2.Totalvariation

TotalVariation(TV)isawidelyusedcostfunctioninimage pro-cessing[14,5,13,17].Bounded TV set was introduced by [3,2] for denoisingproblems:

CT V

= {

x

|

T V

(

x

)

≤



}

(10)

Fig. 1.1and2balls.

whichis thesetofimageswhoseTV isbelowagiven



.Theset CT V isalsoaclosedandconvexset.Thereforeitcanbeusedinany ProjectionontoConvexSets(POCS)baseddeconvolutionproblem.

Inthisarticleweintroducetheboundedfilteredvariation(FV) set,whichisbasedonthefollowingcostfunction:

F V

(

x

)

=



n1,n2

|

x

[

n1

,

n2

] − (

g

∗

x

)

[

n1

,

n2

]|

(11)

whereg isalow-passfilter.Anylow-passfiltercanbeusedinEqn.

(11).In1-D TVfunction g

[

n

]

issimplyequal to

δ

[

n

−

1

]

because 1-D TVfunctionis T V

(

x

)

=



_n

|

x

[

n

]

−

x

[

n

−

1

]|

.The boundedFV setisdefinedasfollows:

CF V

= {

x

|

F V

(

x

)

≤



}

(12)

ItcanbeshownthatCF V isalsoaclosedandconvexset. We perform projections onto the set CF V in an approximate mannerintwosteps.Letxk bethecurrentimagethatwewantto projectontothesetCF V.Theimagexk isdividedintoitslow-pass andhigh-passcomponentsxk,lo

=

xk

∗

g andxk,hi

=

xk

−

xk,lo using thelow-passfilterg,respectively.Weprojectthehigh-passfiltered componentontothe



1-ballandobtain:

xk,hp

=

P1

(

xk,hi

)

(13)

where P1 representtheorthogonalprojection operationontothe



1-ball.Finally,wecombinethelow-passcomponentoftheimage

withxk,hp toobtaintheapproximateprojectionontothesetCF V:

xk+1

=

xk,lo

+

xl,hp (14)

wherex_k+1 containsthelow-passcomponentsoftheoriginal

im-agebutitshigh-passcomponentsareregulatedbythe



1-ball.

1.3. Boundedenergyset

Both



1-balland



2-ballarewellknownsetsusedinimage

re-constructionproblems.The



1-ballis

C1

=





n1,n2

|

x

[

n1

,

n2

]| ≤



(15) andthe



2-ballis

C2

=





n1,n2

|x[n

1

,

n2

]|

2

≤



2

(16) InFig. 1



1and



2ballsin R2areshown.

Orthogonalprojectionontothe



2 ballisimplementedby

sim-plyscalingthecoefficientsoftheinputvector xj.Letthe

||

xj

||

2

=

α

2_and

_α

2

_>

_

2_{.Theprojectionx}

j+1 ofxj ontothe



2ballisgiven

by

xj+1

=



(



2

_/

_α

2

₎

_x

j (17)

Thevectorxjcanbeprojectedontothe



1ballusingDuchietal.’s

method [4]. But it can alsobe approximatelyprojected onto the



1 ballintwo steps. Wefirstproject thecurrentiterateonto the



2 ballandthentothe



1 ball.Forthesecondstep,theprojection

Fig. 2. Flow chart of the proposed algorithm.

vector x_j+1 obtainedbythe aboveequationcanbeprojectedonto

theboundaryhyperplaneofthe



1 ballasfollows:

xj+2

=

xj+1

−



n1,n2 xj+1

[

n1

,

n2

] +



N (18)

Some ofthe entries ofxj+1 mayleave the quadrant or RN after

theaboveprojectionoperation.Weforcethosecoefficientsofx_j+1

tozero.

1.4. Iterativedeconvolutionmethod

The POCS based iterative deconvolution algorithm uses the above projectionoperationsina successivemanner. Letxn bethe n-thiterate.Thenextiteratex_n+1 isobtainedby

xn+1

=

PF VP1P2P1,1P1,2

. . .

PL,MPφxn (19)

where PF V representstheprojectionoperationontoCF V,P2

rep-resentstheprojectionoperationontothe



2ball, Pn1,n2 represents

the projection operation onto the hyperplane Cn1,n2 and Pφ

rep-resentsthe projectionoperation ontothehyperplane Cφ,andthe size of the image is L

×

M. Since all of the sets are closed and convextheiteratesconvergetotheintersection

Cint

=

C1

∩

C2

∩

n1,n2Cn1,n2

∩

Pφ (20)

regardless oftheinitial inputimage x1 providedthat theset Cint is non-empty. The order ofprojections inEq. (19) isimmaterial. The convergence speed ofiterations canbe improvedby linearly combiningtheorthogonalprojections[3,10].

Fig. 3. (a)projDeconvoutputasitisprocessedbyMATLAB2016a.Theinputisblurredwithadiskfilterwithradius=9pixels.TheaddedGaussiannoisehasσ=0.01. (b) CommandlineprintsprovideinformationofthecurrentiterateandthePSNR.

Fig. 4. projDeconvoutputwithzeronoiseBlurryimage.Theinputisblurredwitha diskfilterwithradius=9pixels.Notethateventhoughnoadditionalnoiseis intro-duced,bordercroppingreducesthequalityofFourierTransformphaseextraction.

Westoptheiterationswhen

||x

n+1

−

xn

|| ≤

γ

(21)

where

γ

isapre-specifiednumberorwestoptheiterationsafter afixednumberofsteps.

Theproposedmethodistestedinvariousdatasetsandthe re-sultsarepresentedinSection3.Inthetests,itisobservedthatthe flatbackgroundimageswithlesssignificant Fouriertransform in-formationinhigherfrequenciesshowmeaningfulimprovementsin termsofrecoveredimagequality.Alsotherecoveredimagequality is higher in highly blurry images compared to the other well-knowndeconvolutionmethodsmentionedinthispaper.These re-sultsareparalleltotheintendeduseofthesoftwareinmicroscopic imagingandMPI.

2. Softwaredescription

Presented software DeconvApp v2 is a MATLAB function im-plemented on the 2016a version of the program. MATLAB is an engineeringtool that has versions for Windows, Linux andMAC

Fig. 5. Deconvolutionresultsofanimagewithhightotalvariationvalues.ThePSNR resultofthepresentedsoftware(20.40 dB)isnotasgoodastheresultsofWiener (20.60 dB)orLucy–Richardsondeconvolutions(20.51 dB).Inthisimagethe back-groundcontains“random”objects.Asaresultthephaseoftheimageappearstobe random.Thismaybethecauseofinferiorperformance.

environments. Thesoftwarewritten makes useoftheImage pro-cessing toolbox functions such as conv2 and fft2, therefore this toolboxisrequiredinordertorunthesoftware.

In orderto make aquick visual comparisonwithwell known

Wiener and Lucy–Richardson deconvolution methods, software

alsodisplaysthedeconvolutionresultsofthesealgorithms. A demo script which calls blurs an image and calls the pre-senteddeconvolutionfunctionisalsopresented.Thesoftwarethat generates the results presented in this paperis also available in Github.

2.1. Softwarearchitecture

Theapplicationrunsinaniterativemanner.As inputweneed a blurryimage andthe blurringfilterandalsoa limit forthe it-erations. The iteration number limit is used to break iteration if reached,howeverwealsostopiterationsiftheconsecutiveiterates are similar to each other. The similarity measure is the

mean-squareerrorofresultingimagesineachiterationtobelowerthan agiventhreshold.

Wealsorequiretheoriginalimage asaninput inorderto cal-culateperformancefactorPSNR.

Anoutline ofthealgorithmandrelatedsoftwareisgivenasa flowchartinFig. 2.

Procedurestartsby extractingestimate Fourierphase informa-tionfromblurryimage.Theinitialestimate fortheoriginalimage canbearandomimage.Startingfromthisinitialestimate,Fourier

Fig. 6. projDeconvoutputwithGaussiannoise(σ=0.01)Blurryimage.Theinputis image12003fromBSDS500datasetblurredwithadiskfilterwithradius=13px.

transformphaseisrecoveredfromblurryimageinFourierdomain andthespatialdomainfeatures arerestoredby



1 ball projection

iteratively.

3. Illustrativeexamples

In the GitHub repository[20],a demo script ‘demo.m’is pre-sented. This script blurs a test image ‘barbara.png’ andcalls the deconvolutionfunctionfordeblurringtheresultingimage. Alterna-tively proj Deconv can be calleddirectly foranyblurry grayscale image.Whenwecallthedemoscriptawindowpopsupasshown inFig. 3.a.

Inthecommandwindow, itispossibletotrackthecurrent it-erationindexandthequalityofdeconvolutionintermsofPSNRas showninFig. 3.b.

Forblurringtheinputimage,weusedafunctioncalled ‘blurIm-age’ which is also available in GitHub. We used thisfunction to crop resulting filtered image to prevent border gradient which naturally occurs in convolution operations. With removal of this bordereffectwe bettersimulatereallifeblurryimages.Thismay introduce additionalcomplexities indeconvolutionproblems. The difficultiescreatedby croppededgesare observableintheresults ofwell knowndeconvolutionmethodssuch asWienerandLucy– Richardson. However, our software deals with thisby artificially generatingagradientusingtheinputfilter.

We tested the algorithm with zero noise case as shown in

Fig. 4.Forthiswemanuallymodifiedthecode‘blurImage’andset std_dev variableto0.

Additionally wetestedoursoftwareusingthedatasetpresented in[7]andtheBSDS500dataset[11].Thesedatasetsareusedin de-convolutionandsegmentationapplications[8,16].Thefirstdataset thatwastakenfrom[7]consistsoffourimagesin

.

mat format.The filesincludegrayscaleimagematricesaswellasapsfandblurred image. TheBSDS500 datasetincludes500colorimageswhichwe convertedtograyscaleinordertouseastestdata.Weselected20 images fromthis datasetwhich satisfies a low noise background feature thatisessential formeaningfulFourierphaseinformation forrecovery.Anexampleoflow qualityreconstruction isgivenin

Fig. 5. We tested these selected imagesagainst regularized filter deblurringfunctionofMatlab(‘deconvreg’)[6]also.

Asourzerophaseassumptionrequires,weneededtoapplyour filterstotheoriginal imagesinthedatasets.Weapplied diskand Gaussianfilterswithchangingradius/varianceontotheimages.We alsoappliedaGaussiannoisewith

σ

=

0

.

01.Itisimportanttonote thattheimagepixelvaluerangeis0to1.

Anexamplerunofthesoftwarewithanimagefromdataset[7]

isgiveninFig. 6.

We compared the results of our algorithm with well known Wiener and Lucy–Richardsondeconvolution algorithms. We used MATLABimplementationsofthesealgorithms.ForLucy–Richardson we onlyusedtheblurry image andthefilter,andfortheWiener case,wesetthenoiseparametertobe0.01.Thecomparisonof re-sults of Barbara, Lenaand the dataset of [7] for disk filterwith disksizes6,9and12aregivenintheTable 1.

Table 1

PSNRcomparisonofdeconvolutionalgorithmsfordiskshapefilterwithvariousradiusr.Firstfourimagesare fromdatasetpresentedin[7]andthelasttwoimagesarestandardtestimages.

Presented software (dB) Wiener deconvolution (dB) Lucy–Richardson method (dB)

r: 6px 9px 12px 6px 9px 12px 6px 9px 12px Im1 22.14 20.88 19.76 19.66 17.99 16.72 20.28 18.91 17.75 Im2 22.25 21.48 20.77 21.91 20.47 19.59 22.55 20.85 19.84 Im3 22.04 21.20 20.23 22.53 20.99 19.74 23.43 21.52 20.06 Im4 25.01 23.29 22.35 21.61 19.89 18.64 22.44 20.49 18.97 Barbara 21.67 21.15 20.65 21.40 20.41 19.73 22.28 21.32 20.18 Lena 24.56 23.80 23.04 23.71 22.29 21.14 25.30 23.56 22.37

Table 2

PSNRcomparisonofdeconvolutionalgorithmsforGaussianshapefilter.Firstfourimagesarefromdataset pre-sentedin[7]andthelasttwoimagesaretestimages.

Presented software (dB) Wiener deconvolution (dB) Lucy–Richardson method (dB)

σ: 4 7 11 4 7 11 4 7 11 Im1 21.30 20.05 19.84 21.76 17.87 16.49 22.36 19.15 16.94 Im2 21.44 20.62 20.75 22.26 19.88 19.36 22.37 20.48 19.61 Im3 21.21 20.12 20.12 23.15 20.23 19.54 23.06 20.80 20.07 Im4 23.80 22.24 22.02 23.67 19.58 18.23 24.03 20.60 18.09 Barbara 21.34 20.75 20.70 22.16 20.24 19.52 22.41 21.20 20.53 Lena 23.99 23.10 23.08 25.14 21.94 20.80 25.60 23.41 21.96 Table 3

PSNRcomparisonofdeconvolutionalgorithmsover20imagesfromBSDS500datasetfordiskfilterofradius:13px.

Image nr. in BDSD500 dataset Presented software (dB) Wiener deconvolution (dB) Lucy–Richardson method (dB) Regularized filter deconvolution (dB)

3096 27.98 25.91 27.31 20.21 12003 21.20 17.45 18.28 14.08 15001 21.51 17.25 19.46 12.90 15088 21.72 19.54 20.10 15.29 19021 21.54 21.14 21.01 15.12 22013 21.54 14.21 13.98 13.38 24004 20.04 19.65 19.44 15.57 29030 23.38 19.45 20.42 15.86 35049 24.45 24.64 23.89 18.27 41096 21.86 20.25 20.59 15.39 48017 17.75 17.62 17.47 12.94 69000 20.48 18.09 18.75 13.57 69015 23.78 21.67 22.25 17.08 100007 23.82 17.92 20.10 13.83 107045 21.43 20.08 20.58 16.85 135037 23.09 16.70 19.21 11.11 153077 24.18 21.36 21.97 13.29 226022 21.93 20.94 21.08 18.05 260081 19.61 16.27 17.19 13.93 302003 21.95 20.53 20.27 10.79 Table 4

MeanPSNRcomparisonofdeconvolutionalgorithmsover20imagesfromBSDS500dataset. Disk filter (dB) Gaussian filter (dB)

r: 6px r: 9px r: 12px σ: 4 σ: 7 σ: 11 Proposed software 24.10 22.99 22.16 23.04 21.23 20.03 Wiener deconvolution 21.83 20.49 19.53 22.96 20.67 19.07 Lucy–Richardson method 22.67 21.12 20.17 23.08 21.06 19.63 Regul. filter deconvolution 19.60 17.86 15.65 19.24 17.15 16.38

We performed deconvolution over Gaussian blurred images

fromthe same sets. The results of thesedeconvolution tests are giveninTable 2.

A larger test is applied to selected images from BSDS500 database. In this casewe compared the results of the proposed method with the results of Wiener, Lucy–Richardson and Regu-larizedfilterdeblurring.Forthelastmethod,Matlab function ‘de-convreg’isused.Table 3presentsasmallportionoftheextensive resultsofthetestedimagesfor13px-radiusdisk-typeblurring fil-ter.

Statistical resultsofall the testsmadefor differenttypesand scalesofblurringpsfsaregivenTable 4.

FromTables 1 and2,weobservethatastheamountofblur in-creases,reductioninoutput quality intermsofPSNR istheleast comparedtoothertwowellknowndeconvolutionmethods. There-foresuccessratesseemtosurpassforlargerpsfs.InTable 4wesee thatforflatbackgroundimagestheproposedmethodoutperforms well-knowndeconvolutionmethods.

4. Impact

The software provided with this paper and its various ver-sions(Other algorithmsthat we developedwhich relyonFourier

transformphase)areusedinourgroupfordeblurringmicroscopic images.Recentlywe applieda similarconcept toMPIwhich pro-duces highly blurry images that are unusable without a proper deconvolution.Weobservedthateventhoughwellknown deblur-ringalgorithmssuchasWienerandLucy–Richardsonmethodsmay be used in such cases,for highblurconditions, an addition ofa phase condition shows great improvement. Therefore we believe theprovided softwaremaybe preferable inhighblur deconvolu-tionproblemssuchasMPI.

Theproposed methodandaccompanying softwareaimsto de-convolveblurryimagesusingprojectionsontoconvexsets(POCS). While the



1–



2 ball projectionsandfiltered variation isused in

regularization purposes,the projectiononto Fourierphase is per-formed to correctany phase deformationsduring iterations. The known zerophase psf requirementis a mustforthe software to proceedcorrectly.

However the algorithm can be extended to non-zero phase

cases by subtracting thephase of the blurringfunction fromthe observedimagetoestimatethephaseoftheoriginalimage.

Itisknownthatrapidlychangingimageswithhightotal varia-tion has noise-likebehavior in their Fourier domain counterpart, resulting in high frequency components in Fourier domain and random lookingphase informationafternoise additionand

crop-ping. Insuch casesitmaynot be always feasibletoestimate the phaseoftheoriginalimage inFourierdomain.Thereforeitis rec-ommendedto useoursoftware inordinary imagesincluding mi-crobiologyimagesandMPIimages.

The presented software and its other versions has not been usedincommercialsettings.Theprovidedsoftwareistheonly ver-sion that is madepublicly available. It is madeopen source and thereforeitispossibleto modifythecodeaccordingtotheneeds oftheapplication.

5. Conclusions

Inthispaper,weproposeaPOCSbaseddeconvolutionsoftware thattakesadvantageofFourierTransformPhasecharacteristics of symmetricblurs.We alsouseorthogonalprojectionsonto



1 and



2ballsandFilteredVariationinordertoregularizetherestoration

process in each iteration. In addition, we introduce a computa-tionally efficient method to estimate projection onto the



1 ball

withoutcomplicatedorderingoperations.Orderingtheimage data isan O rder

(

N

)

sortingoperationbut N

=

L2 foran imageofsize L by L.

Ourapproximateprojectionconsistsoftwosteps.Firstwemake anorthogonalprojectionontothe



2 ballwhichisbasicallyscaling

the data by a constant. Then we make an orthogonal projection onto theboundary hyperplane of the



2 ball whichis an almost

multiplicationfreeoperation.

Theexperimentalresultsshowthattheprovidedsoftware out-performswellknownWienerandLucy–Richardsonmethodsinflat backgroundimageswithlargeblurpsfs.

Appendix A. Supplementarymaterial

Supplementarymaterialrelatedtothisarticlecanbefound on-lineathttps://doi.org/10.1016/j.dsp.2017.11.004.

References

[1]PatrizioCampisi,KarenEgiazarian,BlindImageDeconvolution:Theoryand Ap-plications,CRCPress,2016.

[2]GiovanniChierchia,NellyPustelnik,J-C.Pesquet,BéatricePesquet-Popescu, Epi-graphicalprojectionandproximaltoolsfor solvingconstrainedconvex opti-mizationproblems,SignalImageVideoProcess.9 (8)(2015)1737–1749.

[3]PatrickL.Combettes,J.-C.Pesquet,Imagerestorationsubjecttoatotalvariation constraint,IEEETrans.ImageProcess.13 (9)(2004)1213–1222.

[4]JohnDuchi,ShaiShalev-Shwartz,YoramSinger,TusharChandra,Efficient pro-jectionsontothel1-ballforlearninginhighdimensions,in:Proceedingsofthe 25thInternationalConferenceonMachineLearning,ACM,2008,pp. 272–279.

[5]PascalGetreuer,TotalvariationdeconvolutionusingsplitBregman,Image Pro-cess.OnLine2(2012)158–174.

[6] TheMathWorksInc.Deblurringimagesusingaregularizedfilter,2017.[Online; accessed25-August-2017].

[7]AnatLevin,YairWeiss,FredoDurand,WilliamT.Freeman,Understandingand evaluatingblinddeconvolutionalgorithms,in:ComputerVision and Pattern Recognition,2009.CVPR2009.IEEEConferenceon,IEEE,2009,pp. 1964–1971.

[8]AnatLevin,YairWeiss,FredoDurand,WilliamT.Freeman,Understandingblind deconvolutionalgorithms,IEEETrans.PatternAnal.Mach.Intell.33 (12)(2011) 2354–2367.

[9]Yanjun Li,KiryungLee, YoramBresler, Identifiabilityinblinddeconvolution with subspaceorsparsityconstraints,IEEETrans.Inf. Theory62 (7)(2016) 4266–4275.

[10]RobertJ.MarksII,HandbookofFourierAnalysis&ItsApplications,Oxford Uni-versityPress,2009.

[11]D.Martin,C.Fowlkes,D.Tal,J.Malik,Adatabaseofhumansegmentednatural imagesanditsapplicationtoevaluatingsegmentationalgorithmsand measur-ingecologicalstatistics,in:Proc.8thInt’lConf.ComputerVision,vol. 2,July 2001,pp. 416–423.

[12]AlanV.Oppenheim,JaeS.Lim,Theimportanceofphaseinsignals,Proc.IEEE 69 (5)(1981)529–541.

[13]DanielePerrone,PaoloFavaro,Aclearerpictureoftotalvariationblind decon-volution,IEEETrans.PatternAnal.Mach.Intell.38 (6)(2016)1041–1055.

[14]LeonidI.Rudin,StanleyOsher, EmadFatemi,Nonlineartotalvariationbased noiseremovalalgorithms,Phys.D:NonlinearPhenom.60 (1)(1992)259–268.

[15]M.I.Sezan,HenryStark,Imagerestorationbythe methodofconvex projec-tions:part2-applicationsandnumericalresults,IEEETrans.Med.Imaging1 (2) (1982)95–101.

[16]LibinSun,SunghyunCho,JueWang,JamesHays,Edge-basedblurkernel es-timationusingpatchpriors,in:ComputationalPhotography(ICCP),2013IEEE InternationalConferenceon,IEEE,2013,pp. 1–8.

[17]MohammadTofighi,KivancKose,A.EnisCetin,Denoisingimagescorruptedby impulsivenoiseusingprojectionsontotheepigraphsetofthetotalvariation function(pes-tv),SignalImageVideoProcess.9 (1)(2015)41–48.

[18]Mohammad Tofighi, Onur Yorulmaz, Kivanç Köse, Deniz Cansen Yıldırım, RengülÇetin-Atalay,A.EnisCetin,Phaseandtvbasedconvexsetsforblind deconvolutionofmicroscopicimages,IEEEJ.Sel.Top.SignalProcess.10 (1) (2016)81–91.

[19]H.Trussell,M.Civanlar,TheLandweber iterationandprojectionontoconvex sets,IEEETrans.Acoust.SpeechSignalProcess.33 (6)(1985)1632–1634.

[20] OnurYorulmaz,A.EnisCetin,Githubrepositoryforsubmittedarticle,https:// github.com/onuryorulmaz/Deconvolution.git,2016.

[21]OnurYorulmaz,O.BurakDemirel,TolgaCukur,EmineU.Saritas,A.EnisCetin, Ablinddeconvolutiontechniquebasedonprojectionontoconvexsetsfor mag-neticparticleimaging,IEEETrans.Biomed.Eng.(2017),submittedfor publica-tion.

[22]DanC.Youla,HeywoodWebb,Imagerestorationbythemethodofconvex pro-jections:part1:theory,IEEETrans.Med.Imaging1 (2)(1982)81–94.

OnurYorulmaz receivedhisElectricalandElectronicsEngineeringB.Sc. degreefromBilkentUniversityin2009.HereceivedhisM.S.E.degreefrom BilkentUniversityin2012.HisM.S.E.thesisisconcernedwithfoodsafety applicationsusingimageprocessingmethods.HeiscurrentlyaPh.D. can-didateinBilkentUniversity.

His researchinterestsare applications ofdigital imageenhancement methodsonbiomedicalimages.

A.EnisCetin studied ElectricalEngineering atthe OrtaDoguTeknik Universitesi (METU).AftergettinghisB.Sc. degree,hegothisM.S.E. and Ph.D. degrees inSystemsEngineering fromthe MooreSchoolof Electri-calEngineeringattheUniversityofPennsylvaniainPhiladelphiain1987. Between1987–1989,hewas AssistantProfessorofElectricalEngineering attheUniversityofToronto,Canada.SincethenhehasbeenwithBilkent University,Ankara, Turkey.Currently heis afull professor. During sum-mers of 1988, 1991, 1992 he was with Bell Communications Research (Bellcore),NJ,USA.Hespent1994–1995academicyearatKocUniversity inIstanbul,and1996–1997academicyearattheUniversityofMinnesota, Minneapolis,USA asavisitingassociateprofessor.Hewas onsabbatical leaveas avisitingprofessor atUC SanDiegoin2016. Heiscurrentlya research professor inthe Dept. ofElectrical and Computer Engineering, UniversityofIllinoisatChicago,Chicago,Illinois.

Prof.CetinisafellowofIEEE.HewasanAssociateEditoroftheIEEE TransactionsonImageProcessingbetween1999and2003,EURASIP Jour-nal ofAppliedSignal Processing (Springer), andSignal Processing (Else-vier)and amemberoftheSPTMtechnicalcommitteeoftheIEEESignal Processing Society. He founded the Turkish Chapter of the IEEE Signal Processing Society in1991. He was Signal Processing and AES Chapter Coordinator in IEEE Region-8 in 2003. He received the young scientist award ofTUBITAK(Turkish Scientificand Technical ResearchCouncil)in 1993.Hewastheco-chairoftheIEEE-EURASIPNonlinearSignalandImage Processing Workshop(NSIP’99)whichwasheldinJune1999inAntalya, Turkey. Hewas also the technical co-chairof the EuropeanSignal Pro-cessingConference,EUSIPCO-2005.HewasontheeditorialboardsofIEEE SignalProcessingMagazine,IEEECASforVideoTechnology,EURASIP Jour-nals,SignalProcessingandJournalofAdvancesinSignalProcessing(JASP). Currently,heistheEditor-in-ChiefofSignal,Imageand VideoProcessing SIViP,Springer.HeholdsfiveUSpatents.