CN107657605A - A kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint - Google Patents
A kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint Download PDFInfo
- Publication number
- CN107657605A CN107657605A CN201710810756.6A CN201710810756A CN107657605A CN 107657605 A CN107657605 A CN 107657605A CN 201710810756 A CN201710810756 A CN 201710810756A CN 107657605 A CN107657605 A CN 107657605A
- Authority
- CN
- China
- Prior art keywords
- mrow
- sieve plate
- msub
- pixel
- bmo
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
- G06T7/0012—Biomedical image inspection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/50—Depth or shape recovery
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10101—Optical tomography; Optical coherence tomography [OCT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30041—Eye; Retina; Ophthalmic
Landscapes
- Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Health & Medical Sciences (AREA)
- General Health & Medical Sciences (AREA)
- Medical Informatics (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Radiology & Medical Imaging (AREA)
- Quality & Reliability (AREA)
- Image Analysis (AREA)
Abstract
The invention discloses a kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, this method is using surface segmentation method before the Bruch films opening point assay method based on k mean clusters and active profile and the sieve plate based on energy constraint, initial profile first using the dendrogram of k mean clusters as C V active contour models, extract the Bruch films opening point in profile, the area-of-interest of surface segmentation before sieve plate is obtained further according to the position of opening point, the preceding surface of sieve plate is partitioned into using energy constraint method, case depth before sieve plate is finally gone out according to the structure measurement of two steps.This method acquired results are better than existing method, it is and consistent with expert's craft calibration result, it can solve the problem that expert needs the problem of wasting time and energy of case depth before demarcation measurement sieve plate manually in clinical diagnosis, there is positive effect to the early screening and clinical diagnosis of glaucoma.
Description
Technical field
The invention belongs to image identification technical field, is related to a kind of preceding surface of the sieve plate based on active profile and energy constraint
Depth measurement method.
Background technology
Glaucoma is second-biggest-in-the-world blinding disease, and it to the peripapillary optic ganglion cell aixs cylinder of people by breaking
It is bad, cause the visual field of people to lack.Due to the irreversibility of glaucoma, morning detection, early discovery and early treatment to glaucoma can
Slow down the process of disease development.But because the pathogenesis of glaucoma is not yet completely clearly, therefore to the risk factors of glaucoma
Research be still current hot issue.
Lamina Cribrosa (lamina cribrosa) is located at sclera, is that retinal ganglion cell axon passes eyes
Position, in sieve texture, it is responsible for providing nutrition and structural support to retinal ganglion cell axon, and help to maintain intraocular
And outside eye between barometric gradient.The anatomical structure of sieve plate is as follows:Ten thousand retinal ganglion cell axons of about 1.2-1.5
It is gathered in and regards nipple, eyes is partly left by the inside (i.e. Bruch's films opening) and outside (i.e. sclera) of canalis opticus.Depending on
Nerve is a part for central nervous system, and outside is surrounded by three layers of envelope (endocranium, the spider web continued by three layers of meninx
Film, pia mater), endocranium cavity of resorption and cavum subarachnoidale are also extended to around optic nerve therewith, and cavum subarachnoidale is full of cerebrospinal fluid.
At sclera metapore, 1/3 forms thin layer fiber every there is many sieve apertures to allow retinal ganglial cellses in dividing plate before scleral tissue
Aixs cylinder is walked, and rear 2/3 connects with endocranium.Research shows the sieve plate risk factors maximum with glaucoma --- intraocular pressure has one
Determine relation, would generally be along with the rise of intraocular pressure in the pathogenic process of glaucoma, and the rise of intraocular pressure can cause simultaneously
The rearward displacement of sieve plate and flexural deformation etc..Wang Ning profits et al. propose optic ganglion cell axonal injury and across sieve plate pressure differential
(difference of intraocular pressure and intracerebral pressure) is relevant, and sieve plate can also show the feature of rearward displacement under this pressure differential.
Clinic of means of optical coherence tomography (Optic Coherence Tomography, the OCT) technology in ophthalmology
Using only 20 years, its technological innovation was rapid, now had become one of most important inspection of clinical ophthalmology.It is right that the technology passes through
Tissue emissions coherent light, reclaims the reflected light and scattering light of tissue, and is combined with the delay of time, can get tissue
The information of two-dimentional fault information or overall retina three-dimensional image.Except in real time monitoring, it is noninvasive the advantages that outside, OCT is most important
Feature is its high resolution, and the fine structure observed is cross-sectional structure, meets pathological normal observation custom, is
The live body retinal morphology research of ophthalmology researcher provides technical support.The OCT course of work is as follows:Sent by low-coherence light source
Low coherence light, two-beam line is divided into by interferometer, it is a branch of to enter detection light path, direct incident intraocular, by intraocular difference group
The interface knitted fires back, to provide the thickness and range information of the various tissues of intraocular;Another beam enters with reference to light path, by
Know that the reference mirror of space length reflects.Two-beam is integrated into one in optical fiber coupler, due to being transmitted into coming for reference mirror
The distance returned distance and be transmitted into the given structure of intraocular accurately matches, therefore produces interference, is detected by light-sensitive detector.Adjust
Signal input computer after solution carries out computing, obtains the optical coherence tomography image of testee.Because eye inner tissue has
Different depth and space structure, so meeting generation time is poor between two-beam line, this time difference is called the optical delay time.
The time difference that detector is calculated is related to tolerance principle using low-coherent light, the information of Tissue reflectance can be obtained.Get
After these information, one-dimensional scan image information is calculated by computer, usually the row in two dimensional image.
In recent years, the progress of OCT technology gradually solves the problems, such as that deep tissues can not be imaged.Heidelberg company OCT's
Strengthen Depth Imaging (Enhanced Depth Imaging, EDI) pattern can by light collection in the deeper position on eyeground,
Such as sclera, choroid and sieve plate.This imaging technique have beneficial to the pathologic of live body sieve plate change such as sieve plate to
Backward shift etc. is further studied, and improves the accuracy of glaucoma early diagnosis.We are usually using case depth before sieve plate
To represent the rearward displacement of sieve plate, the measuring method of this sieve plate depth is:With Bruch films opening point (Bruch ' s
Membrane Opening, BMO) composition BMO reference planes (sieved with the preceding surface of sieve plate as reference plane, measurement datum
The distance between interface between being organized before plate tissue and sieve plate).Therefore, the measurement of case depth is divided into two steps before sieve plate:
BMO is determined and surface segmentation before sieve plate.After obtaining the position on BMO reference planes and the preceding surface of sieve plate, case depth can before sieve plate
To be drawn by measuring the distance between two faces.
Because the research to sieve plate just starts to walk, some other existing methods can not also be advantageously applied to sieve plate figure
In handling, so related achievement in research is relatively fewer.
A) BMO is determined
BMO points are one of most important biomarkers in OCT image, and it represents optic disk edge in OCT faultage images
Position, therefore it is widely used in optic disk segmentation, disk is the problems such as measurement along.
2013, Fu et al. proposed a kind of BMO assay methods rebuild based on dictionary learning and low-rank matrix, this method
Be partitioned into first retina internal limiting membrane (ILM) layer and retinal pigment epithelium (RPE) layer and divide the image into training region and
Candidate region, low-rank dictionary then is trained with the training region of image, and goes out whole image with the data reconstruction of candidate region,
Obtain reconstruction error curve, including luminance errors, local binary patterns (LBP) error and apart from biasing.Finally by three mistakes
Poor curve draws the position of the position i.e. BMO points of optic disk.The shortcomings that this method, is the BMO that can not effectively determine glaucoma lesion
Position.
2014, Akram et al. is proposed a kind of determined BMO position based on the model of deconvolution.The model is thought
Each layer in OCT image can be modeled by two compositions, be respectively:With a curve come model the skeleton of layer and with one or
One group of wave filter models the thickness of layer.As long as the parameter of the two compositions is obtained, it is possible to uniquely determine every on OCT image
One layer of position.The model is entered using Monte Carlo Markov Chain (Monte Carlo Markov Chain, MCMC) method
Row parametric measurement, reach higher accuracy rate, but also bring the problem of less efficient simultaneously.
2015, Wang et al. proposed a kind of BMO assay method, the method use five kinds of methods and determines respectively
BMO position, shape facility of the first method based on ILM layers, texture of remaining four kinds of method based on BMO neighbouring positions are special
Sign.After obtaining the result of five kinds of methods, select wherein three immediate values and ask their average value to be surveyed as final BMO
Determine result.The problems such as less efficient, robustness is poor be present in this method.
Hussian et al. proposes a kind of BMO assay methods based on graph theory, and this method measures three reference layers first
Position, it is retinal nerve fiber (RNFL) layer, external plexiform layer (ONL) layer and RPE layers respectively.Then the position using layer and spy
The weight as figure is levied, BMO position is drawn using a kind of graph search algorithm, has reached higher accuracy rate.
B) surface segmentation before sieve plate
The preceding surface of sieve plate is sieve plate tissue and the interface organized before sieve plate, not solid because the sieve plate of people comes in every shape
Fixed shape facility, contrast unobvious sometimes can also be traditional based on gradient because the disease such as glaucoma causes to lack
Feature, the dividing method effect of color and shape and bad, some existing Hierarchical Segmentation methods are not suitable for table before sieve plate yet
In the segmentation of face.
2015, Akram et al. proposed the automatic division method on the preceding surface of the first sieve plate, and this method is based on Ma Erke
Husband's random field and MCMC sampling algorithms, classification problem is converted into by segmentation problem, it is only necessary to knows whether each point belongs to sieve
The preceding surface of plate.In units of each point, extract the feature of the point and its neighbours' point establish a Markov with
Airport, the parameter of markov random file is then determined using MCMC algorithms, and sample by iteration to obtain final surface.
However, above-mentioned existing method generally exist it is less efficient, to glaucoma lesion image not robust the problems such as, cause it
Practicality in the examination and clinical diagnosis of glaucoma it is poor.
The content of the invention
The invention provides a kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, its purpose
It is, overcomes the problem of efficiency is low in the prior art and robustness is not strong.
A kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, comprises the following steps:
Step 1:After gray level image to be detected is carried out into inverse processing, all pixels point in image is handled to inverse and adopted
Clustered with k means clustering methods, obtain dendrogram;
Step 2:The energy function of C-V active contour models is constructed, C-V active profile dies are used as by the use of the profile of dendrogram
The initial curve of type, using the curve during energy function minimums of C-V active contour models as objective contour, with objective contour
The terminating point of lower boundary obtains BMO and refers to projection line as BMO points, two BMO points of connection;
Step 3:Handle and chosen in image positioned at BMO with reference to the region below projection line from reflection, it is corresponding with selected areas
Area-of-interest of the minimum rectangular area as surface segmentation before sieve plate, wherein, the coboundary of selected areas is BMO with reference to throwing
The projection line in the horizontal direction of hachure;
Step 4:The preceding surface candidate point of one sieve plate, and foundation are extracted from area-of-interest based on energy function in each row
Integrity constraint, rejecting processing is carried out to all preceding surface candidate points of sieve plate, obtain controlling point set;Control point set is entered
Row curve matching obtains floor lever of sieve tray contour line;
Step 5:Projection line and floor lever of sieve tray contour line are referred to based on BMO, calculate case depth before sieve plate;
Using BMO with reference to the midpoint of projection line and its nasal side and two points of each 100 microns of temporo side as measurement point, at three
Vertical line is made with reference to projection line to BMO in measurement point and intersected with surface profile line before sieve plate, the average value of three length of perpendicular is
Case depth before sieve plate to be measured.
Further, the dendrogram is M (i, j, k), is to handle inverse to belong to the minimum cluster of pixel value in image
The pixel of the class of center representative retains, and the gray value of rest of pixels point is set to 0:
Wherein, cmAnd ckM-th and k-th of cluster centre, and c are represented respectivelymRepresent the current cluster of pixel (i, j)
Center, I (i, j) be image in pixel (i, j) pixel value, I (cm) and I (ck) it is respectively cluster centre cmAnd ckPlace
Pixel value.
Further, the energy function of the C-V active contour models is:
E (C)=μ L (C)+v*Area (inside (C))+λ ∫inside(C)(I(x,y)-Ii)2dxdy
+ω∫outside(C)(I(x,y)-Io)2dxdy
Wherein, I (x, y) is the pixel value of the pixel (x, y) in dendrogram, and L (C) represents length of a curve, Area
(inside (C)) represents the area that contour curve is surrounded, IiAnd IoThe respectively pixel of outer and in profile the pixel of profile
Value;Inside (C) represents the pixel in profile, and outside (C) represents the pixel outside profile in dendrogram;μ joins for length
Number, span is (0,1);V is area parameters, value 0;λ and ω is respectively the weight coefficient of self-energy and outer energy, λ
=ω=1.
Further, according to following during the energy function based on surface segmentation before sieve plate respectively arranges from area-of-interest
The preceding surface candidate point of one sieve plate of formulas Extraction:
Wherein, SjThe preceding surface candidate point of sieve plate in jth row is represented, n represents total columns of pixel in area-of-interest;Ij
(i) the brightness value of the ith pixel point of jth row pixel in area-of-interest, gradient (I are representedj(i) I) is representedj(i)
Gradient Features value, α and β represent the brightness of pixel and the weight parameter of Gradient Features respectively, and α span is
[0.3,0.4], β span is [0.5,0.6], and arg max () represent to find the position for the point for making energy function maximum.
Total columns of pixel is the width of region of interest area image in area-of-interest;
Arg max () are to cause α * Ij(i)+β*gradient(Ij(i) when) obtaining maximum, corresponding pixel is made
For Sj;
Further, the integrity constraint is:
Wherein, Control (j) represents control point set, and d represents offset distance control parameter, and span is (10,15);S
For whole candidate point ordinate set, YjThe ordinate of the candidate point on j-th of A-scan is represented, average (S) waits to be whole
The average value of reconnaissance set ordinate.
Further, using B-spline fitting algorithm point set will be controlled to carry out curve fitting to obtain floor lever of sieve tray contour line.
Beneficial effect
The present invention proposes a kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, this method
Based on active contour model, extract the BMO points in OCT image and BMOBMO refers to projection line, reuse energy constraint method
It is partitioned into the preceding surface of sieve plate.Finally using the result of surface segmentation before BMO measure and sieve plate, case depth before sieve plate is measured.
This method has used the active contour model based on region in BMO points measure, avoids making measure using traditional Gradient Features
As a result by the abnormal influence such as blood vessel shade, lesion;Before sieve plate in surface segmentation, this method employs the think of of energy function
Want solve the problem of preceding surface of sieve plate is without the fixed traditional characteristic such as shape, size, and sieve is fitted using curve matching
The lack part of plate;Test result indicates that the accuracy rate peace of this method surface segmentation before the BMO accuracy determined and sieve plate
Existing method, the result more demarcated manually close to expert are superior in slip.The present invention can solve the problem that expert in clinical diagnosis
When need manually demarcation measurement sieve plate before the problem of wasting time and energy of case depth, early screening and clinical diagnosis to glaucoma
With positive effect.
Brief description of the drawings
Fig. 1 is the flow chart of the present invention, wherein, (a) is its logical flow chart, and (b) is that it performs the flow chart of step;
Fig. 2 is the schematic diagram of Measurement of surface deepth before sieve plate;
Fig. 3 is the BMO measurement results of this method, wherein, (a) is measurement result Fig. 1 of this method;(b) to be corresponding with (a)
The result demarcated by hand of expert;(c) it is measurement result Fig. 2 of this method;(d) expert corresponding with (c) demarcates by hand
As a result;(e) it is measurement result Fig. 3 of this method;(f) result demarcated by hand for expert corresponding with (e);
Fig. 4 is surface segmentation result before the sieve plate of this method;
The contrast schematic diagram that Fig. 5 is demarcated by hand for case depth before the sieve plate of this method measurement with expert, wherein, (a) is
Measurement result Fig. 1 of this method;(b) result demarcated by hand for expert corresponding with (a);(c) it is the measurement result of this method
Fig. 2;(d) result demarcated by hand for expert corresponding with (c);(e) it is measurement result Fig. 3 of this method;(f) it is right with (e)
The result that the expert answered demarcates by hand.
Embodiment
Below in conjunction with accompanying drawing, the present invention will be further described.
The present invention proposes a kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, including four
Individual step:First with initial profile of the dendrogram of k mean algorithms as C-V active contour models, then contouring lower boundary
Terminating point, as BMO points, then extract in image comprising BMO with reference to projection line region below minimum rectangle as sieve
The area-of-interest of surface segmentation before plate, candidate point is extracted according to energy function to each row of interesting image regions, will
The candidate point for not meeting constraints is left out, and is partitioned into whole surface with curve-fitting method, finally calculates table before sieve plate
Face depth.
Shown in idiographic flow such as Fig. 1 (a), shown in implementation procedure such as Fig. 1 (b) of this method, its specific implementation step is as follows:
Step 1:After gray level image to be detected is carried out into inverse processing, all pixels point in image is handled to inverse and adopted
Clustered with k means clustering methods, obtain dendrogram;
Image to be detected I (m, n) after given inverse, k mean algorithms are divided into m × n pixel in image predefined
K class, it is smaller according to the otherness of pixel in each class class, and the principle that otherness between class and class is larger, Wo Menke
The problem of making cost function (representing the otherness between pixel) minimum so that problem is converted into, then obtain every a kind of cluster
Center so that the otherness of each point and nearest cluster centre is minimum.Used cost function is as follows:
Wherein c is predefined classification number, d (i, j, k)=| | I (i, j)-I (ck)||2For two points on description image it
Between diversity factor metric function because OCT image has level feature, therefore use the pixel value tag of image as difference
Property measurement, I (i, j) be image in pixel, ckFor cluster centre.So make the minimum dendrogram M (i, j, k) of cost function
It is calculated according to lower formula.
Wherein, cmAnd ckM-th and k-th of cluster centre, and c are represented respectivelymRepresent the current cluster of pixel (i, j)
Center, I (i, j) be image in pixel (i, j) pixel value, I (cm) and I (ck) it is respectively cluster centre cmAnd ckPlace
Pixel value.The pixel for belonging to the class that the minimum cluster centre of pixel value represents in figure is retained, remaining sets to 0, as to be measured
Determine dendrogram.After obtaining dendrogram, the central point per a kind of all pixels is exactly new cluster centre.
Step 2:The energy function of C-V active contour models is constructed, C-V active profile dies are used as by the use of the profile of dendrogram
The initial curve of type, using the curve during energy function minimums of C-V active contour models as objective contour, with objective contour
The terminating point of lower boundary obtains BMO and refers to projection line as BMO points, two BMO points of connection;
Because the result of active contour model is largely influenceed by given initial curve, therefore by upper one
Step obtained by k mean clusters figure as initial curve C, the C-V active contour model of active contour model energy function such as
Following formula:
E (C)=μ L (C)+v*Area (inside (C))+Einside+Eoutside
=μ L (C)+v*Area (inside (C))+λ ∫inside(C)(I(x,y)-Ii)2dxdy
+ω∫outside(C)(I(x,y)-Io)2dxdy
Wherein, I (x, y) is the pixel value of the pixel (x, y) in dendrogram, and L (C) represents length of a curve, Area
(inside (C)) represents the area that contour curve is surrounded, IiAnd IoThe respectively pixel of outer and in profile the pixel of profile
Value;Inside (C) represents the pixel in profile, and outside (C) represents the pixel outside profile in dendrogram;μ joins for length
Number, its value is determined by the size of target in image, if the size of target is larger, μ value is also larger, otherwise μ
Value will be smaller, span 0-1, generally takes μ=0.5.V is area parameters, and λ and ω are respectively self-energy and outer energy
Weight coefficient.Generally take λ=ω=1, v=0.C-V models obtain the distribution of its energy by constructing energy function, with energy
Function ever-reduced method promotes curve to approach the edge of target, and most target is split from background at last.It breaks away from
Dependence in classical edge partitioning algorithm such as LOG and traditional Snake models to image gradient, have to OCT image good
Segmentation ability.
Step 3:Handle and chosen in image positioned at BMO with reference to the region below projection line from reflection, it is corresponding with selected areas
Area-of-interest of the minimum rectangular area as surface segmentation before sieve plate, wherein, the coboundary of selected areas is BMO with reference to throwing
The projection line in the horizontal direction of hachure;
After obtaining positions of the BMO with reference to projection line, area-of-interest is used as with reference to the floor projection line of projection line using BMO
The width of coboundary, i.e. rectangular area;Region lower boundary is image base, and right boundary is respectively the position of two BMO points.
Step 4:The preceding surface candidate point of one sieve plate, and foundation are extracted from area-of-interest based on energy function in each row
Integrity constraint, rejecting processing is carried out to all preceding surface candidate points of sieve plate, obtain controlling point set;Control point set is entered
Row curve matching obtains floor lever of sieve tray contour line;
Brightness and Gradient Features around the preceding surface of sieve plate build the energy function of surface segmentation before sieve plate, and
The preceding surface candidate point of a sieve plate, such as following formula are selected according to energy function in each A-scan (i.e. each row, similarly hereinafter):
Wherein, SjThe preceding surface candidate point of sieve plate in jth row is represented, n represents total columns of pixel in area-of-interest;Ij
(i) the brightness value of the ith pixel point of jth row pixel in area-of-interest, gradient (I are representedj(i) I) is representedj(i)
Gradient Features value, α and β represent the brightness of pixel and the weight parameter of Gradient Features respectively, and α span is
[0.3,0.4], β span is [0.5,0.6], and arg max () represent to find the position for the point for making energy function maximum
Put.After obtaining candidate point set, because the morphological feature of sieve plate varies with each individual, and the influence of blood vessel shade, there are some times
Reconnaissance is in fact, it is necessary to find out and deleted the common trait of this part candidate point not on the position on the preceding surface of sieve plate
Go.Integrity constraint is as follows:
Wherein, Control (j) represents control point set, and d represents offset distance control parameter, is needed with this state modulator
The excessive candidate point of the deviant of removal, value is between 10-15.S is whole candidate point ordinate set, YjRepresent j-th
The ordinate of candidate point on A-scan, average (S) are the average value of whole candidate point set ordinate.In order to ensure
Smooth surface before the sieve plate arrived, we are carried out curve fitting on obtained control point set Control (j), and the fitting used is calculated
Method is B-spline fitting algorithm.
Step 5:Projection line and floor lever of sieve tray contour line are referred to based on BMO, calculate case depth before sieve plate;
Obtain before BMO measurement results and sieve plate after surface segmentation result, with BMO with reference to the midpoint of projection line and its nasal side and
Two points of each 100 microns of temporo side as measurement point, in three measurement points to BMO with reference to projection line make vertical line and with before sieve plate
Intersect on surface.Such as Fig. 2, the average value of three length of perpendicular is case depth before sieve plate to be measured.
To verify the performance of method proposed by the present invention, 30 width images of 18 samples are gathered as data set, and with special
Ground Truth figures corresponding to the preceding surface conduct of BMO points and sieve plate that family demarcates by hand.Based on described in step 3, α and β parameters
Whether effective method has been largely fixed it, has been found through experiments that, when α increases, curve can be close to sieve plate region;When β increases
When big, curve can be close to internal limiting membrane region.For the robustness of ensuring method, α=0.3-0.4, β=0.5-0.6 are set.
The BMO assay methods and NLSC sieve plates that the present invention proposes the computational methods of proposition and Hussian in 2015 et al.
Preceding surface segmentation method is compared.In BMO measure, we are with the mean error and variance of measuring point and Ground Truth
As evaluation criterion;Before sieve plate in surface segmentation, we are then using error rate as evaluation criterion.
Table 1BMO measurement results
From table 1 it follows that method proposed by the present invention (48.17 microns of mean error) in accuracy is better than
The BMO assay methods (54.18 microns of mean error) that Hussian et al. is proposed, and in stability, the standard deviation of this method
It is 51.32 microns, is also slightly better than its 53.74 microns.Fig. 3 lists this method and obtains segmentation result and expert's manual segmentation
BMO points and BMO reference planes, as can be seen from Figure 3 due to avoiding Gradient Features, the present invention can be very good to avoid by big
The problem of image blurring that blood vessel shadow band comes, can also obtain accurately splitting very much on the EDI-OCT images more than blood vessel shade
As a result, it is and consistent with the result that expert demarcates manually.
(error be present with Ground Truth in the error rate before sieve plate in surface segmentation using the preceding surface sampling point of sieve plate
Point the ratio between number and total sample number) evaluation criterion is used as, error rate grade is defined as from 0 to 2:It is small that grade 0 represents mean error
The ratio between total sample number is accounted in the point of 3 pixels, grade 1 represents point of the mean error between 3 pixels and 5 pixels and accounts for sample
The ratio between sum, grade 2 represent mean error and account for the ratio between total sample number more than the point of 5 pixels.By point of proposition method of the present invention
Cut result and a kind of BFPS edges dividing method and the preceding surface of currently the only sieve plate based on Markov chain Monte-Carlo
Automatic division method --- NLSC methods are contrasted, as a result as shown in table 2.
Surface segmentation result before the sieve plate of table 2
100-150 evaluation point is up-sampled in every pictures.As can be seen from Table 2 compared with Ground Truth,
In error rate grade 0 and grade 1, method proposed by the invention respectively reaches 76.5% and 17.6%, higher than NLSC methods
73.7% and 16.1% and the 40.9% of BFPS methods and 22.8%, and in the higher error rate grade 2 of error, institute of the present invention
Moving party's rule only has 5.9%, is better than the 36.3% of the 10.2% and BFPS methods of NLSC methods.Institute of the present invention as can be seen here
The method of proposition is better than NLSC methods and BFPS methods on the whole.Fig. 4 lists pair of experimental result and Ground Truth
Than, it can be seen that the result that experimental result is demarcated manually with expert is very close on image, and curve is also very smooth, from
Also there is reasonability from the point of view of in histology.
Sieve plate depth measured by last the inventive method is compared with Ground Truth, still using mean error
It is as shown in table 3 as evaluation criterion, its result with standard deviation.
The preceding surface automatic measurement result of the sieve plate of table 3
It can be seen that, the sieve plate depth and the resultant error of expert's demarcation that the present invention measures are smaller from table 3, basic one
Causing, Fig. 5 gives more intuitively result, it can also be seen that the craft demarcation of this method and expert are basically identical from figure, and
Accuracy rate is higher, and it is 22.1 microns averagely to have error in label, and no error in label is 28.7 microns.
Specific embodiment described herein is only to spirit explanation for example of the invention.Technology belonging to the present invention is led
The technical staff in domain can be made various modifications or supplement to described specific embodiment or be replaced using similar mode
Generation, but without departing from the spiritual of the present invention or surmount scope defined in appended claims.
Claims (6)
1. a kind of Measurement of surface deepth method before sieve plate based on active profile and energy constraint, it is characterised in that including following
Step:
Step 1:After gray level image to be detected is carried out into inverse processing, all pixels point in image is handled to inverse and uses k
Means clustering method is clustered, and obtains dendrogram;
Step 2:The energy function of C-V active contour models is constructed, C-V active contour models are used as by the use of the profile of dendrogram
Initial curve, using the curve during energy function minimums of C-V active contour models as objective contour, with the following of objective contour
The terminating point on boundary obtains BMO and refers to projection line as BMO points, two BMO points of connection;
Step 3:Handle and chosen in image positioned at BMO with reference to the region below projection line from reflection, with corresponding to selected areas most
Area-of-interest of the small rectangular area as surface segmentation before sieve plate, wherein, the coboundary of selected areas refers to projection line for BMO
Projection line in the horizontal direction;
Step 4:The preceding surface candidate point of one sieve plate is extracted in each row from area-of-interest based on energy function, and according to complete
Property constraints, rejecting processing is carried out to all preceding surface candidate points of sieve plate, obtain control point set;Point set march will be controlled
Line is fitted to obtain floor lever of sieve tray contour line;
Step 5:Projection line and floor lever of sieve tray contour line are referred to based on BMO, calculate case depth before sieve plate;
Using BMO with reference to the midpoint of projection line and its nasal side and two points of each 100 microns of temporo side as measurement point, in three measurements
Vertical line is made with reference to projection line to BMO on point and intersected with surface profile line before sieve plate, the average value of three length of perpendicular is to be measured
Case depth before the sieve plate of amount.
2. according to the method for claim 1, it is characterised in that the dendrogram is M (i, j, k), is to scheme inverse processing
The pixel for belonging to the class that the minimum cluster centre of pixel value represents as in retains, and the gray value of rest of pixels point is set to 0:
<mrow>
<mi>M</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mi>I</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mi>f</mi>
<mo>|</mo>
<mo>|</mo>
<mi>I</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>I</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>c</mi>
<mi>m</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msub>
<mo>|</mo>
<mn>2</mn>
</msub>
<mo>&le;</mo>
<mo>|</mo>
<mo>|</mo>
<mi>I</mi>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>I</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>c</mi>
<mi>k</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msub>
<mo>|</mo>
<mn>2</mn>
</msub>
<mo>,</mo>
<mi>m</mi>
<mo>&NotEqual;</mo>
<mi>k</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mn>0</mn>
</mtd>
<mtd>
<mrow>
<mi>o</mi>
<mi>t</mi>
<mi>h</mi>
<mi>e</mi>
<mi>r</mi>
<mi>s</mi>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
</mrow>
Wherein, cmAnd ckM-th and k-th of cluster centre, and c are represented respectivelymThe current cluster centre of pixel (i, j) is represented,
I (i, j) be image in pixel (i, j) pixel value, I (cm) and I (ck) it is respectively cluster centre cmAnd ckThe pixel at place
Value.
3. according to the method for claim 2, it is characterised in that the energy function of the C-V active contour models is:
E (C)=μ L (C)+v*Area (inside (C))+λ ∫inside(C)(I(x,y)-Ii)2dxdy+ω∫outside(C)(I(x,y)-
Io)2dxdy
Wherein, I (x, y) is the pixel value of the pixel (x, y) in dendrogram, and L (C) represents length of a curve, Area (inside
(C) area that contour curve is surrounded, I) are representediAnd IoThe respectively pixel value of outer and in profile the pixel of profile;inside
(C) pixel in profile is represented, outside (C) represents the pixel outside profile in dendrogram;μ is length parameter, value model
Enclose for (0,1);V is area parameters, value 0;λ and ω is respectively the weight coefficient of self-energy and outer energy, λ=ω=1.
4. according to the method for claim 1, it is characterised in that the energy function based on surface segmentation before sieve plate is from sense
In interest region the preceding surface candidate point of a sieve plate is extracted in each row according to below equation:
<mrow>
<msub>
<mi>S</mi>
<mi>j</mi>
</msub>
<mo>=</mo>
<mi>arg</mi>
<munder>
<mrow>
<mi>m</mi>
<mi>a</mi>
<mi>x</mi>
</mrow>
<mi>i</mi>
</munder>
<mrow>
<mo>(</mo>
<mi>&alpha;</mi>
<mo>*</mo>
<msub>
<mi>I</mi>
<mi>j</mi>
</msub>
<mo>(</mo>
<mi>i</mi>
<mo>)</mo>
<mo>+</mo>
<mi>&beta;</mi>
<mo>*</mo>
<mi>g</mi>
<mi>r</mi>
<mi>a</mi>
<mi>d</mi>
<mi>i</mi>
<mi>e</mi>
<mi>n</mi>
<mi>t</mi>
<mo>(</mo>
<mrow>
<msub>
<mi>I</mi>
<mi>j</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>)</mo>
</mrow>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
<mo>...</mo>
<mi>n</mi>
</mrow>
Wherein, SjThe preceding surface candidate point of sieve plate in jth row is represented, n represents total columns of pixel in area-of-interest;Ij(i) table
Show the brightness value of the ith pixel point of jth row pixel in area-of-interest, gradient (Ij(i) I) is representedj(i) ladder
Spend characteristic value, α and β represent the brightness of pixel and the weight parameter of Gradient Features respectively, α span for [0.3,
0.4], β span is [0.5,0.6], and argmax () represents to find the position for the point for making energy function maximum.
5. according to the method for claim 4, it is characterised in that the integrity constraint is:
<mrow>
<mi>C</mi>
<mi>o</mi>
<mi>n</mi>
<mi>t</mi>
<mi>r</mi>
<mi>o</mi>
<mi>l</mi>
<mrow>
<mo>(</mo>
<mi>j</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>S</mi>
<mi>e</mi>
<mi>t</mi>
<mo>{</mo>
<mi>j</mi>
<mo>|</mo>
<mfrac>
<mrow>
<msub>
<mi>dY</mi>
<mi>j</mi>
</msub>
</mrow>
<mrow>
<mi>d</mi>
<mi>j</mi>
</mrow>
</mfrac>
<mo>*</mo>
<mfrac>
<mrow>
<msub>
<mi>dY</mi>
<mrow>
<mi>j</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
<mrow>
<mi>d</mi>
<mi>j</mi>
</mrow>
</mfrac>
<mo>&le;</mo>
<mn>0</mn>
<mo>,</mo>
<mo>|</mo>
<msub>
<mi>Y</mi>
<mi>j</mi>
</msub>
<mo>-</mo>
<mi>a</mi>
<mi>v</mi>
<mi>e</mi>
<mi>r</mi>
<mi>a</mi>
<mi>g</mi>
<mi>e</mi>
<mrow>
<mo>(</mo>
<mi>S</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<mo><</mo>
<mi>d</mi>
<mo>}</mo>
<mo>,</mo>
<mi>j</mi>
<mo>=</mo>
<mn>1</mn>
<mo>...</mo>
<mi>n</mi>
</mrow>
Wherein, Control (j) represents control point set, and d represents offset distance control parameter, and span is (10,15);S is whole
Individual candidate point ordinate set, YjThe ordinate of the candidate point on j-th of A-scan is represented, average (S) is whole candidate point
Gather the average value of ordinate.
6. according to the method described in claim any one of 1-5, it is characterised in that will control point set using B-spline fitting algorithm
Carry out curve fitting to obtain floor lever of sieve tray contour line.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710810756.6A CN107657605B (en) | 2017-09-11 | 2017-09-11 | A kind of sieve plate front surface depth measurement method based on active profile and energy constraint |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710810756.6A CN107657605B (en) | 2017-09-11 | 2017-09-11 | A kind of sieve plate front surface depth measurement method based on active profile and energy constraint |
Publications (2)
Publication Number | Publication Date |
---|---|
CN107657605A true CN107657605A (en) | 2018-02-02 |
CN107657605B CN107657605B (en) | 2019-12-03 |
Family
ID=61128062
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710810756.6A Active CN107657605B (en) | 2017-09-11 | 2017-09-11 | A kind of sieve plate front surface depth measurement method based on active profile and energy constraint |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN107657605B (en) |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109528155A (en) * | 2018-11-19 | 2019-03-29 | 复旦大学附属眼耳鼻喉科医院 | A kind of intelligent screening system and its method for building up suitable for the concurrent open-angle glaucoma of high myopia |
CN110288540A (en) * | 2019-06-04 | 2019-09-27 | 东南大学 | A kind of online imaging standards method of carbon-fibre wire radioscopic image |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103236065A (en) * | 2013-05-09 | 2013-08-07 | 中南大学 | Biochip analysis method based on active contour model and cell neural network |
EP2967317A2 (en) * | 2013-03-14 | 2016-01-20 | Carl Zeiss Meditec AG | Multimodal integration of ocular data acquisition and analysis |
CN106373168A (en) * | 2016-11-24 | 2017-02-01 | 北京三体高创科技有限公司 | Medical image based segmentation and 3D reconstruction method and 3D printing system |
-
2017
- 2017-09-11 CN CN201710810756.6A patent/CN107657605B/en active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
EP2967317A2 (en) * | 2013-03-14 | 2016-01-20 | Carl Zeiss Meditec AG | Multimodal integration of ocular data acquisition and analysis |
CN103236065A (en) * | 2013-05-09 | 2013-08-07 | 中南大学 | Biochip analysis method based on active contour model and cell neural network |
CN106373168A (en) * | 2016-11-24 | 2017-02-01 | 北京三体高创科技有限公司 | Medical image based segmentation and 3D reconstruction method and 3D printing system |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109528155A (en) * | 2018-11-19 | 2019-03-29 | 复旦大学附属眼耳鼻喉科医院 | A kind of intelligent screening system and its method for building up suitable for the concurrent open-angle glaucoma of high myopia |
CN110288540A (en) * | 2019-06-04 | 2019-09-27 | 东南大学 | A kind of online imaging standards method of carbon-fibre wire radioscopic image |
CN110288540B (en) * | 2019-06-04 | 2021-07-06 | 东南大学 | Carbon fiber wire X-ray image online imaging standardization method |
Also Published As
Publication number | Publication date |
---|---|
CN107657605B (en) | 2019-12-03 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Bogunović et al. | RETOUCH: The retinal OCT fluid detection and segmentation benchmark and challenge | |
Kafieh et al. | A review of algorithms for segmentation of optical coherence tomography from retina | |
US10194866B2 (en) | Methods and apparatus for reducing artifacts in OCT angiography using machine learning techniques | |
Baroni et al. | Towards quantitative analysis of retinal features in optical coherence tomography | |
Sigal et al. | Eye-specific IOP-induced displacements and deformations of human lamina cribrosa | |
WO2022007352A1 (en) | Three-dimensional choroidal vessel imaging and quantitative analysis method and apparatus based on optical coherence tomography system | |
Norman et al. | Dimensions of the human sclera: thickness measurement and regional changes with axial length | |
Esmaeelpour et al. | Choroidal Haller's and Sattler's layer thickness measurement using 3-dimensional 1060-nm optical coherence tomography | |
US7992999B2 (en) | Automated assessment of optic nerve head with spectral domain optical coherence tomography | |
Knott et al. | Spatial correlation of mouse photoreceptor-RPE thickness between SD-OCT and histology | |
Xu et al. | Automated geographic atrophy segmentation for SD-OCT images based on two-stage learning model | |
Guo et al. | A framework for classification and segmentation of branch retinal artery occlusion in SD-OCT | |
Zhang et al. | A novel technique for robust and fast segmentation of corneal layer interfaces based on spectral-domain optical coherence tomography imaging | |
Hussain et al. | An automated method for choroidal thickness measurement from Enhanced Depth Imaging Optical Coherence Tomography images | |
Dhaini et al. | Automated detection and measurement of corneal haze and demarcation line in spectral-domain optical coherence tomography images | |
CN107657605B (en) | A kind of sieve plate front surface depth measurement method based on active profile and energy constraint | |
CN109744996B (en) | OCT image BMO position positioning method | |
Eladawi et al. | Optical coherence tomography: A review | |
Monemian et al. | Mathematical analysis of texture indicators for the segmentation of optical coherence tomography images | |
EP3417401B1 (en) | Method for reducing artifacts in oct using machine learning techniques | |
Piorkowski et al. | Selected aspects of corneal endothelial segmentation quality | |
US20180315185A1 (en) | System and method to analyze various retinal layers | |
Mari et al. | A digital staining algorithm for optical coherence tomography images of the optic nerve head | |
Cheung et al. | Region of interest densitometry analysis of Descemet membrane endothelial keratoplasty dehiscence on anterior segment optical coherence tomography | |
Ibrahim et al. | Volumetric quantification of choroid and Haller's sublayer using OCT scans: An accurate and unified approach based on stratified smoothing |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |