CN106373165B - Tomography composograph method for reconstructing and system - Google Patents
Tomography composograph method for reconstructing and system Download PDFInfo
- Publication number
- CN106373165B CN106373165B CN201610796748.6A CN201610796748A CN106373165B CN 106373165 B CN106373165 B CN 106373165B CN 201610796748 A CN201610796748 A CN 201610796748A CN 106373165 B CN106373165 B CN 106373165B
- Authority
- CN
- China
- Prior art keywords
- mrow
- voxel
- view
- value
- msub
- 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.)
- Active
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T11/00—2D [Two Dimensional] image generation
Abstract
The present invention relates to a kind of tomography composograph method for reconstructing and system.The above method includes：S10, the non-update area for obtaining the n-th 1 3-D views, calculate the voxel difference of m reference pictures and the n-th 1 3-D view, interpolation fitting is carried out to non-update area, the compensation base value S20 of non-each voxel of update area is obtained, obtains the n-th 1 compensation images；S30, the data for projection to the n-th angle are rebuild, and obtain the n-th 3-D view；S40, judge whether n is equal to n0, if it is not, n then is updated into n+1, return to step S10, if so, then entering step S50；S50, using the n-th 3-D view as m reconstruction image results, judge whether m is equal to default iterations, if it is not, into step S60, if so, then using m reconstruction images result as final tomography composograph；S60, processing is weighted to m reference pictures according to m reconstruction images result and iterations m, obtains m+1 reference pictures, and m is updated to m+1, return to step S10.
Description
Technical field
The present invention relates to tomography composograph processing technology field, more particularly to a kind of tomography composograph method for reconstructing
And system.
Background technology
Mammary gland X ray photographing system is the important means of diagnosis of breast disease and examination, and it utilizes plate for forcing and supporting plate
The measured targets such as breast are pressed, and the irradiated tissue by grenz ray to above-mentioned measured target, by detector to wearing
Saturating X ray is received and handled, the two-dimensional projection image of measured target structure after being oppressed.However, above-mentioned mammary gland X is penetrated
In line camera chain, the three-dimensional structure of measured target can not avoid the overlapping and mutual coverage of tissue on the two-dimensional projection image, and
And it can not still obtain the structural information to longitudinal depth of the measured targets such as breast.
To realize the acquisition of above-mentioned measured target three-dimensional information and solving the problems such as tissue therein is overlapping, grind at present
Send the mammary gland digital tomosynthesis system (Digital based on digital tomography synthetic technology (Digital Tomosynthesis, DTS)
Breast Tomosynthesis, DBT), such system sets up x-ray bulb on the basis of the mammary gland X ray photographing system of routine
Independent rotational movement structure, in the case where compression device, measured target (breast etc.), detector keep rigid condition, by bulb along circle
Arc orbital rotational motion, multiple low dose exposure is carried out to the tissue in the measured target of pressing in limited angular range, so
Several detector cells on detector are read afterwards and obtain a series of two-dimensional projection datas, are finally generated using image rebuilding method
The measured target 3-D view of pressing.The 3-D view that DBT is obtained can solve two-dimensional digital mammography image (Digital
Mammogram, DM) tissue is overlapping and the problems such as false positive, and can be realized with the DM images under equal oppression state good
Contrast and combination, the diagnosis efficiency of mammary gland disease is improved to a certain extent.
If directly carrying out DBT image reconstructions using analytic reconstruction methods such as filtered back projection FBP and FDK, lack sampling has
It is more difficult to wave filter design to limit sparse angular projection, the missing of projection information is difficult to obtain the tomography composite diagram being accurately satisfied with
Picture.And traditional backprojection reconstruction (BP) algorithm, though basic three-dimensional structure information can be obtained quickly, reconstructed results lack essence
Thin details, clinic can not be directly applied to.
Parallel algebraic reconstruction technique (Simultaneous Algebraic Reconstruction Technique,
SART it is) that current DBT rebuilds relatively accurate effective method.However, in DBT imaging systems, the x-ray bulb of circular motion and
The angle that detector normal is formed is bigger, and longitudinal depth information of acquisition is abundanter.But angle crosses conference causes image objects
The subregion in region (Field of Vision, FOV) exceeds the effective range of detector, so as to cause projection information to block
And missing, projection absent region is formed, and the size of the projection absent region can be with the change of scanning angle and dimension of object
And change.When being updated operation to 3-D view voxel under all angles successively in SART process of reconstruction, corresponding angle
The voxel of lower projection absent region will be unable to obtain back projection's renewal, form non-update area, and the projection of adjacent angular is blocked
Difference will cause each non-update area adjacent edges pixel value in corresponding tomography composograph to occur repeatedly drastically being mutated, this
Kind phenomenon will ultimately result in the tomography composograph obtained in reconstruction and a plurality of banding artifact occurs along scanning direction lateral side regions, from
And the quality and effect that the tomography composograph for influenceing measured target is rebuild.
The content of the invention
Based on this, it is necessary to the technical problem of measured target tomography composograph quality is influenceed for traditional reconstruction scheme,
A kind of new tomography composograph method for reconstructing and system are provided.
A kind of tomography composograph method for reconstructing, comprises the following steps：
S10, the non-update area for obtaining the (n-1)th 3-D view, calculate m reference pictures and (n-1)th 3-D view
Voxel difference, interpolation fitting is carried out to non-update area according to voxel difference, obtains the benefit of non-each voxel of update area
Repay base value；Wherein, n is more than 1 and is less than or equal to n0Integer, n0To launch the child's hair twisted in a knot-childhood number of degrees of X ray, m is current iteration time
Number, initial value are set as 1；
S20, according to it is described compensation base value compensation, root are modified to the non-update area voxel value of the (n-1)th 3-D view
The voxel value obtained according to correction-compensation obtains the (n-1)th compensation image；
S30, using parallel algebraic reconstruction technique according to (n-1)th compensation image the data for projection of the n-th angle is rebuild,
Obtain the n-th 3-D view；
S40, judge whether n is equal to n0, if it is not, n then is updated into n+1, return to step S10, if so, then entering step
S50；
S50, using the n-th 3-D view as m reconstruction image results, judge whether m is equal to default iterations, if
It is no, into step S60, if so, then using m reconstruction images result as final tomography composograph；
S60, processing is weighted to m reference pictures according to m reconstruction images result and iterations m, obtains m+1
Reference picture, and m is updated to m+1, return to step S10.
A kind of tomography composograph reconstructing system, including：
Acquisition module, for obtaining the non-update area of the (n-1)th 3-D view, calculate m reference pictures and described (n-1)th
The voxel difference of 3-D view, interpolation fitting is carried out to non-update area according to voxel difference, it is each to obtain the non-update area
The compensation base value of voxel；Wherein, n is the integer more than 1 and less than or equal to n0, and n0 is the child's hair twisted in a knot-childhood number of degrees of transmitting X ray, and m is
Current iteration number, initial value are set as 1；
Correcting module, for being repaiied according to the compensation base value to the non-update area voxel value of the (n-1)th 3-D view
Positive compensation, the voxel value obtained according to correction-compensation obtain the (n-1)th compensation image；
Module is rebuild, for the data for projection using parallel algebraic reconstruction technique according to the (n-1)th compensation image to the n-th angle
Rebuild, obtain the n-th 3-D view；
Judge module, for judging whether n is equal to n0, if it is not, n then is updated into n+1, returns and perform acquisition module, if
It is then to perform update module；
Update module, for using the n-th 3-D view as m reconstruction image results, judging whether m is equal to default iteration
Number, if it is not, setup module is performed, if so, then using m reconstruction images result as final tomography composograph；
Setup module, for m reconstruction images result and iterations m to be weighted into processing to m reference pictures, obtain
M+1 is updated to m+1 reference pictures, and by m, returns and performs acquisition module.
Above-mentioned tomography composograph method for reconstructing and system, the non-update area of the (n-1)th 3-D view can be obtained, calculated
The voxel difference of m reference pictures and the (n-1)th 3-D view, interpolation fitting is carried out to non-update area according to voxel difference, obtained
To the compensation base value of non-each voxel of update area, the voxel value of the (n-1)th 3-D view is carried out according to above-mentioned compensation base value
Correction-compensation, the voxel value obtained according to correction-compensation obtains the (n-1)th compensation image, to be carried out to the data for projection of the n-th angle
Rebuild, obtain the n-th 3-D view, and circulate said process, until obtain the n-th 0 3-D views, regard the n-th 0 3-D views as the
M reconstruction image results, and processing is weighted with this result and obtains m+1 reference pictures, continue to enter the data for projection of all angles
Row iteration is rebuild after iterations reaches preset times, using m reconstruction images result as final tomography composograph；
It can eliminate a plurality of banding artifact of the tomography composograph along scanning direction lateral side regions, improve measured target tomography composite diagram
The reconstruction effect of picture.
Brief description of the drawings
Fig. 1 is the tomography composograph method for reconstructing flow chart of one embodiment；
Fig. 2 is the DBT system schematics of one embodiment；
Fig. 3 is the X ray directive measured target schematic diagram of one embodiment；
Fig. 4 is the digital breast model schematic diagram of one embodiment；
Fig. 5 is the first reference picture schematic diagram of one embodiment；
The DBT tomography composograph schematic diagrames that Fig. 6 is rebuild by the present invention of one embodiment；
The DBT tomography composograph schematic diagrames that Fig. 7 is rebuild by the traditional scheme of one embodiment；
Fig. 8 is the tomography composograph reconstructing system structural representation of one embodiment.
Embodiment
The embodiment of the tomography composograph method for reconstructing to the present invention and system is made detailed below in conjunction with the accompanying drawings
Description.
With reference to figure 1, Fig. 1 show the tomography composograph method for reconstructing flow chart of one embodiment, comprises the following steps：
S10, the non-update area for obtaining the (n-1)th 3-D view, calculate m reference pictures and (n-1)th 3-D view
Voxel difference, interpolation fitting is carried out to non-update area according to voxel difference, obtains the benefit of non-each voxel of update area
Repay base value；Wherein, n is the integer more than 1 and less than or equal to n0, and n0 is the child's hair twisted in a knot-childhood number of degrees of transmitting X ray, i.e. DBT systems are swept
The total projection number retouched；N initial value is that 2, m is current iteration number, and initial value is set as 1；
DBT systems (mammary gland digital tomosynthesis system) can be with as shown in Fig. 2 wherein plate for forcing, measured target and detector be protected
To hold and be relatively fixed, X ray bulb launches X ray in n0 angle set in advance to measured target successively along arc track,
To be scanned to corresponding measured target, gather several low dose exposure images and carry out tomography composograph reconstruction.Such as Fig. 2
Shown, with detector when some angles (such as larger angle angle) expose, measured target can be projected to outside detector X ray
Portion, cause the projection information in the measured target region to block and lack, form corresponding projection absent region, cause corresponding reconstructed
The loss of learning of image.The n-th three obtained using parallel algebraic reconstruction technique to the projection renewal amendment in the case where carrying out the n-th angle
Image is tieed up, the image has corresponding absent region.Above-mentioned absent region obtains under the angle in the makeover process of data for projection
Less than the renewal amendment of data for projection, corresponding non-update area is formed, because the n-th 3-D view is according to the (n-1)th compensation figure
Data for projection under picture and the n-th angle is updated amendment gained, therefore the non-update area of the n-th 3-D view and the n-th angle
Under data for projection block that absent region is corresponding, i.e. the non-update area of the n-th 3-D view is the projection under the n-th angle
Data block absent region.If the data for projection received to detector directly carries out SART reconstructions, successively under all angles
When being updated operation to 3-D view voxel, the voxel of corresponding absent region will be unable to be updated, and adjacent angular is scarce
Losing area differentiation will cause absent region edge pixel values to occur repeatedly drastically being mutated, and cause in reconstruction image along scanning direction side
There is a plurality of banding artifact in face region.Above-mentioned measured target can include the tested organ of DBT systems scanning.
Fig. 3 show the X ray directive measured target schematic diagram under the (n-1)th angle, if as shown in figure 3, X ray bulb from
The first end of detector obtains the X ray 101 near the end of detector second with being put down where detector 103 to the second end motion
Face angulation 104, (including the first end of measured target 102 is away from spy for the position relationship of measured target 102 and detector 103
The distance of the first end of device 103 is surveyed, including the second end of measured target 102 distance away from the end of detector 103 second), it can also obtain
The size (such as width, height) of measured target 102 is taken, projecting relation and related objective geometric parameter according to light can calculate
The spatial information for blocking absent region under (n-1)th angle, including absent region 105 in figure (do not project to the area of detector
Domain) corresponding to layer scope, line range and row scope of the voxel in measured target, so as to can determine that in the (n-1)th angle reconstruction institute
The non-update area of the (n-1)th 3-D view obtained.
Above-mentioned steps can obtain the calculation of the basic three-dimensional structure of measured target using back projection (BP) method for reconstructing etc.
Method is rebuild to data for projection, to obtain the initial reference image of measured target, i.e. the first reference picture.Utilize back projection
(BP) method for reconstructing can be calculated using the global information of data for projection well, obtain preferable global three-dimensional information and not
Banding artifact occurs, the data for projection under all angles is rebuild using above-mentioned back projection (BP) method for reconstructing, can be with
The basic three-dimensional structure of scanned object (measured target) is obtained, the amendment of non-update area during carrying out SART reconstructions is mended
Repay.
S20, according to it is described compensation base value compensation, root are modified to the non-update area voxel value of the (n-1)th 3-D view
The voxel value obtained according to correction-compensation obtains the (n-1)th compensation image；
Compensation base value can be superimposed in above-mentioned steps on the voxel value of the non-update area of the (n-1)th 3-D view, with reality
The correction-compensation of the non-update area voxel of existing (n-1)th 3-D view.
S30, using parallel algebraic reconstruction technique according to (n-1)th compensation image the data for projection of the n-th angle is rebuild,
Obtain the n-th 3-D view；
In above-mentioned steps, can first compensating (n-1)th image, angularly n carries out forward projection, obtains the n-th estimated projection number
According to, then the processing such as weighted differences is carried out to the data for projection under above-mentioned n-th estimated projection data and the n-th angle, with to preceding to throwing
After the result that shadow and weighted differences obtain carries out back projection, amendment is weighted to (n-1)th compensation each voxel of image, updated
The n-th 3-D view.When data for projection under above-mentioned n-th angle includes the x-ray bombardment measured target of n-th angle, detection
The data for projection that device receives.
S40, judge whether n is equal to n0, if it is not, n then is updated into n+1, return to step S10, if so, then entering step
S50；
S50, using the n-th 3-D view as m reconstruction image results, judge whether m is equal to default iterations, if
It is no, into step S60, if so, then using m reconstruction images result as final tomography composograph；
S60, processing is weighted to m reference pictures according to m reconstruction images result and iterations m, obtains m+1
Reference picture, and m is updated to m+1, return to step S10.
In the present embodiment, the projection angle that exposure can be originated from X ray terminates the projection angle of exposure to X ray, according to
The secondary reconstruction operation carried out to data for projection corresponding to measured target under all angles as described in step S10 to S30, when not up to
In stopping criterion for iteration (generally default iterations, can also use such as other end conditions of image gradient norm) repetition
Process of reconstruction is stated, so that the tomography composograph finally given, keeps higher with the three-dimensional image information of corresponding measured target
Uniformity.Above-mentioned iterations can according to the reconstruction precision of tomography composograph requirement be configured, such as be arranged to 5 or
10 is equivalent.
, can be directly right in first time iterative process for the data for projection of first angle in above-mentioned iterative process
The data for projection of above-mentioned first angle is rebuild, and to obtain the first 3-D view, relation and scan geometry are projected using light
Parameter determines the non-update area of the first 3-D view, is calculated further according to above-mentioned first 3-D view and corresponding non-update area
The first compensation image is obtained, and the weight of the second 3-D view is carried out according to the data for projection under the compensation image and second angle
Build.If in follow-up each secondary iterative process (in addition to for the first time), m reconstruction images that can be obtained by above an iteration process
As a result initial three-dimensional image is used as, and is weighted with m reconstruction images result and m reference pictures by iterations,
M+1 reference pictures are obtained, the data for projection according to above-mentioned initial three-dimensional image successively all angles is carried out such as step S10 extremely
Reconstruction described in S30；M reconstruction images result i.e. obtained by above an iteration process as the first 3-D view, according to
Above-mentioned first 3-D view and the m+1 reference pictures obtained carry out the reconstruction to the second 3-D view.
The tomography composograph method for reconstructing that the present embodiment proposes, it can calculate and obtain not updating for the (n-1)th 3-D view
Region, the voxel difference of m reference pictures and the (n-1)th 3-D view is calculated, non-update area is inserted according to voxel difference
Value fitting, obtains the compensation base value of non-each voxel of update area, according to above-mentioned compensation base value to the (n-1)th 3-D view not more
Each voxel value of new region is modified compensation, and the voxel value obtained according to correction-compensation obtains the (n-1)th compensation image, according to n-th
The data for projection of angle is rebuild, and obtains the n-th 3-D view, and circulates said process, until the n-th 0 3-D views are obtained, will
The n-th 0 3-D views are as m reconstruction image results, and this time result weighting processing obtains m+1 reference pictures, continues to each
The data for projection of individual angle is iterated reconstruction after iteration reaches default iterations, and m reconstruction images result is made
For final tomography composograph；It can eliminate a plurality of banding artifact of the tomography composograph along scanning direction lateral side regions,
Improve the picture quality of measured target tomography composograph and rebuild effect.
In one embodiment, can also include before above-mentioned steps S10：
Backprojection reconstruction is carried out to multiple data for projection that detector receives, the 3-D view that backprojection reconstruction is obtained is set
It is set to the first reference picture；
The data for projection of first angle is rebuild using parallel algebraic reconstruction technique, obtains the first 3-D view；
The non-update area of the first 3-D view is obtained, calculates the voxel of the first reference picture and first 3-D view
Difference, interpolation fitting is carried out to non-update area according to voxel difference, obtains the compensation base value of non-each voxel of update area；
The compensation base value and the first reference picture are modified compensation to the voxel value of the first 3-D view, according to repairing
Just compensate obtained voxel value and obtain the first compensation image；
The data for projection of second angle is rebuild according to the first compensation image using parallel algebraic reconstruction technique, obtained
Second 3-D view.
The non-update area of above-mentioned first 3-D view is the region that detector does not receive data for projection under first angle,
Absent region is blocked in projection i.e. under the angle.The voxel difference of first reference picture and first 3-D view is above-mentioned the
The difference of the voxel value of each voxel voxel value of corresponding voxel respectively and in the first 3-D view in one reference picture.
The present embodiment carries out regional compensation according to the first 3-D view and the first reference picture, and to the projection of second angle
The reconstruction of the second 3-D view is carried out, the accuracy of the second 3-D view of reconstruction can be improved, and is eliminated under first angle not
The voxel mutation of update area.
In one embodiment, above-mentioned multiple data for projection received to detector carry out backprojection reconstruction, by back projection
The process that the 3-D view that reconstruction obtains is arranged to the first reference picture can include：
The data for projection of multiple different angles received to detector is according to detector position, ray source position, object chi
The three-dimensional relationship such as very little and light projection principle, carry out backprojection reconstruction, by back projection to each voxel of target three-dimensional image
Rebuild obtained 3-D view and be arranged to the first reference picture.
When above-mentioned DBT is to measured target transmitting X ray, detector can receive X ray penetrator under different angle
Multiple data for projection of body.
As one embodiment, data for projection is counter under the above-mentioned different angle that the X ray is received to detector is thrown
Shadow is rebuild, and the process that the 3-D view that backprojection reconstruction is obtained is arranged to the first reference picture can include：
The data for projection that the X ray is received to detector carries out the processing of back projection (BP) method for reconstructing, according in difference
The spatial informations such as x-ray source, scanning object (measured target), detector under angle, determine that each voxel X under the angle is penetrated
Line penetrates path, and obtains the pixel coordinate of the projected image under respective path, is sequentially overlapped the angled respective pixel coordinate of institute
The pixel value at place and after weighting processing, you can obtain the basic three-dimensional structure of measured target, the basic three-dimensional structure is set
It is set to the first reference picture.
Above-mentioned back projection (BP) method for reconstructing can be calculated using the global information of data for projection well, be obtained preferably
Global three-dimensional information, be not in banding artifact, the projection received using above-mentioned back projection (BP) method for reconstructing to detector
Data are rebuild, and can obtain the basic three-dimensional structure of scanned object (measured target), as the first reference picture, then
Iterations constantly weights renewal, to carry out the correction-compensation of the non-update area of SART reconstructions.
In one embodiment, above-mentioned steps S30 can include：
To the (n-1)th compensation image, angularly n carries out forward projection, obtains the n-th estimated projection data；
The weighted differences of pixel are carried out according to the data for projection under the n-th estimated projection data and the n-th angle；
After the result back projection processing that above-mentioned forward projection and weighted differences are obtained, to (n-1)th compensation each voxel of image
Amendment is weighted, the n-th 3-D view updated.
The present embodiment can angularly n carries out forward projection to the (n-1)th compensation image, and bulb is sent out under simulation current angular n
Penetrate X ray and penetrate on scanning project objects to detector, calculate each pixel of projected image (detector member) corresponding X ray
Penetrate path, and will penetrate at this successively and n-1 compensation image three-dimensional be corresponded on path spatially each voxel is weighed by respective path
It is weighted again cumulative, obtains the n-th estimated projection data；Entered according to the data for projection under the n-th estimated projection data and the n-th angle
Row weighted differences, the pixel value that estimated projection and acquired projections are corresponded to each pixel carries out difference processing, and presses the pixel
The total length that X ray where (detector member) penetrates scanning object is weighted processing to the pixel difference；To forward projection and
The result that weighted differences obtain carries out back projection's weighting amendment, and bulb, detector, the space for scanning object are angularly corresponded under n
Position, X ray corresponding to projection each coordinate of difference result is obtained, obtain and be pierced corresponding to scanning object as the x-ray path
3 d space coordinate, the path weight value by above-mentioned projection difference result by each voxel on the path of place, weighted superposition to (n-1)th
The pixel value in image corresponding three-dimensional coordinate pixel is compensated, obtains the n-th 3-D view of renewal amendment.
The present embodiment can be mended under the conditions of the DBT that detector finite size and wide-angle scan is rebuild by overall situation projection
The image for repaying amendment carries out SART reconstructions, and reconstruction image scanning direction lateral side regions banding caused by overcoming missing projection data is pseudo-
Shadow, ensure that longitudinal depth information of the tomography composograph of reconstruction, while improve the quality of image reconstruction.
In one embodiment, above-mentioned steps S60 can include：
The voxel value of voxel value and m reference pictures to the m reconstruction image results substitutes into weighting more new formula, meter
Calculate the voxel value of m+1 reference pictures；Wherein, the weighting more new formula is：
Refm+1(i, j, k)=(Refm(i,j,k)+(m-1)pResm(i,j,k))/(1+(m-1)p),
Wherein, Refm+1(i, j, k) represents the voxel value of the i-th row jth row kth layer in m+1 reference pictures, Refm(i,j,
K) voxel value of the i-th row jth row kth layer in m reference pictures, Res are representedm(i, j, k) is represented the in m reconstruction image results
The voxel value of i row jth row kth layers, m represent current iteration number, and p is renewal weight coefficient.Above-mentioned renewal weight coefficient p can be with
It is set greater than 0 value.
The present embodiment is added to m reference pictures by iterations in an iterative process, according to m reconstruction images result
Power renewal, while ensureing the basic structure of initial reference image, the quality of reference picture is stepped up according to iterations, is entered
One step improves the accuracy of the non-update area compensating approach of SART process of reconstruction.In one embodiment, m reference pictures are calculated
Include with the process of the voxel difference of (n-1)th 3-D view：
The voxel value of the voxel value of the m reference pictures and the (n-1)th 3-D view is substituted into voxel differences formula meter respectively
Calculate the voxel difference of the (n-1)th 3-D view；Wherein, the voxel differences formula is：
Diffn-1(i, j, k)=Imgn-1(i,j,k)-Refm(i, j, k),
Wherein, Imgn-1(i, j, k) represents the voxel value of the i-th row jth row kth layer in the (n-1)th 3-D view, Refm(i,j,
K) voxel value of the i-th row jth row kth layer in m reference pictures, Diff are representedn-1(i, j, k) is represented the in the (n-1)th 3-D view
The voxel difference of i row jth row kth layers, m represent current iteration number.
The present embodiment calculates the voxel difference of m reference pictures and the (n-1)th 3-D view, to be inserted to non-update area
Value, the local voxel that can eliminate m reference pictures and the n-th 3-D view of homolographic projection data reconstruction are mutated the mistake brought
Difference.
It is above-mentioned that interpolation fitting is carried out to non-update area according to voxel difference as one embodiment, obtain it is described not more
The process of the compensation base value of each voxel of new region includes：
The voxel difference is substituted into the neighborhood correction value of the non-update area of neighboring mean value formula calculating；
The neighborhood correction value of the non-update area is substituted into the compensation base value of the non-update area of interpolation fitting formula calculating；
The neighboring mean value formula is：
The interpolation fitting formula is：
Wherein, Ndiffn-1(i, k) represents the neighborhood correction value of the i-th row kth layer in the (n-1)th 3-D view, and J represents (n-1)th
The border columns of voxel, JB in 3-D viewn-1(i, k) represents the border of non-update area in (n-1)th 3-D view the i-th row kth
Columns residing for layer, N represent the voxel columns of default neighborhood, Comn-1(i, j, k) represents the i-th row jth in the (n-1)th 3-D view
The compensation base value of row kth layer.
Above-mentioned neighborhood can carry out pre- according to the voxel columns of measured target and the voxel columns of corresponding non-update area
If, it is generally the case that when the data for projection that the n-th angle is carried out to the (n-1)th 3-D view updates, corresponding neighborhood could be arranged to
In (n-1)th 3-D view, from the upper non-update area of an Angles Projections (correspondence that i.e. data for projection of the (n-1)th angle can not update
Region) boundary line start, to the several columns voxel of non-non- update area side；If above-mentioned neighborhood is arranged to never update area side
Boundary line starts, to two row voxels of non-non- update area side, then N=2.
The present embodiment successively can be worked as according to neighborhood columns N, neighborhood correction value Ndiff (i, k) by scanning course bearing fitting
The compensated curve of preceding each voxel of non-update area.Above-mentioned fit approach may include linear fit, it is possible to use the side such as curve matching
Formula calculates the offset of non-each voxel of update area.
In one embodiment, it is above-mentioned that benefit is modified to the voxel value of the (n-1)th 3-D view according to the compensation base value
The process repaid includes：
The voxel value of the compensation base value and the (n-1)th 3-D view is substituted into correction-compensation formula respectively and carries out computing；Institute
Stating correction-compensation formula is：
Wherein, Fixn-1(i, j, k) represents the voxel value of (n-1)th compensation image the i-th row jth row kth layer, JBn-1(i, k) table
Show the border of non-update area in the columns residing for (n-1)th 3-D view the i-th row kth layer, Imgn-1(i, j, k) represents the (n-1)th three
Tie up the voxel value of the i-th row jth row kth layer voxel in image, Comn-1(i, j, k) represents the compensation base of the i-th row jth row kth layer
Value, Refm(i, j, k) represents the voxel value of m reference picture the i-th row jth row kth layers, and m represents current iteration number.
The present embodiment is by the way of superposition is corrected, it is ensured that the non-update area voxel obtained by after superposition is worth to entirely
The compensation of office's projection information, avoid data for projection under traditional scheme and block voxel value mutation caused by missing.
As one embodiment, the above-mentioned voxel by compensation base value superposition 3-D view corresponding with the (n-1)th angle
The process that value is overlapped can also use other stacked systems such as mean filter stacked system.Above-mentioned mean filter superposition can
For：
Wherein, Fixn-1(i, j, k) is the (n-1)th voxel value for compensating image the i-th row jth row kth layer, Com obtained by superpositionn-1
(i, j, k) represents the compensation base value of (n-1)th 3-D view the i-th row jth row kth layer, Refm(i, o, k) represents m reference pictures
The voxel value of i-th row o row kth layers, the length of L mean filters, m represent current iteration number.
Above-mentioned local interpolation method, with reference to the first reference picture using the method for neighborhood fitting and local interpolation to not updating
Region carries out Interpolation compensation, and obtaining information using the global information of data for projection more completely compensates 3-D view, above-mentioned benefit
Repay 3-D view and can be used for the projection renewal operation of SART reconstructions next time, can solve to lack neighborhood border pixel values well
Caused by mutation the problem of gray value acute variation, the generation of banding artifact is avoided.
Tomography composograph method for reconstructing provided by the invention, pass through predetermined angle range scans object to be reconstructed, difference
Corresponding data for projection under each angle is obtained, backprojection reconstruction method processing is carried out to the data for projection of acquisition, obtained scanned
The basic three-dimensional structure of object, and as initial reference image (the first reference picture)；According to mechanical scanning geological information and light
Line projects relation, and data for projection back projection under each angle is handled respectively, accurate point for calculating corresponding non-update area
Cloth；SART reconstructions are carried out to data for projection, under current angular with reference to the first reference picture and to last time angle reconstruction image not
Update area carries out local interpolation, and data correction is carried out to the voxel of non-update area；Carried out according to revised 3-D view
Forward projection weighted differences and back projection's weighting amendment, complete the renewal of 3-D view under Current projection, then complete all angles
Data for projection renewal under degree, is obtained when time iteration result；Reference picture is updated according to iterations and iteration result, continued
Iterative approximation obtains final DBT tomography composographs.Its ensure detector size finite sum wide-angle scanning on the premise of,
Reconstruction image caused by missing projection data can be overcome ensure that DBT images longitudinal direction along scanning direction lateral side regions banding artifact
Depth information, while improve the quality of image reconstruction.
With reference to a concrete application scene, the above-mentioned tomography composograph method for reconstructing of the present invention is illustrated.
As measured target, (model matrix size is 400*600*50 to digital breast model, and voxel is big as shown in Figure 4 for selection
Small is 0.5*0.5*1mm).DBT, which rebuilds (reconstruction of tomography composograph) step, to be included：
1) DBT scannings are carried out to breast model, X ray bulb range of movement is exposed by -30 ° to+30 ° every 3 °
And projection data acquisitions, 21 data for projection Proj are gathered altogethern, n ∈ [1,21], it is assumed that and detector matrix size is 600*800,
Pixel size is 0.4*0.4mm, and corresponding obtained data for projection matrix is 600*800*21；
2) using back projection's (BP) method for reconstructing simultaneously to 21 data for projection Proj of collectionnCarry out backprojection reconstruction,
Reference picture as shown in Figure 5 is obtained, what back projection (BP) method for reconstructing was rebuild to obtain obtains basic three-dimensional structure, above-mentioned to make
First reference image R ef of measured target1；
3) initial three-dimensional image Img is preset0For 0, according to the data for projection Proj of the 1st angle1Carry out SART forward projection
After weighted differences, the 1st reconstruction image Img is obtained after back projection to reconstruction domain weighting amendment1；
4) since the 2nd angle, i.e. during n ＞ 1, principle is projected according to the perspective geometry of setting and light, calculate (n-1)th
3-D view Imgn-1Non- update area, calculate m reference image Rs efm-1With the (n-1)th 3-D view Imgn-1Voxel
Difference, according to voxel difference to non-update area Unn-1Interpolation fitting is carried out, obtains the compensation of non-each voxel of update area
Base value Comn-1；Wherein, n is the integer more than 1 and less than or equal to n0, and n0 is the child's hair twisted in a knot-childhood number of degrees of transmitting X ray, herein n0=
21, m be current iteration number, and initial value is set as 1；
The neighboring mean value formula is：
The interpolation fitting formula is：
Wherein, Ndiffn-1(i, k) represents the neighborhood correction value of (n-1)th 3-D view the i-th row kth layer, and J represents the (n-1)th three
Tie up the border columns of voxel in image, JBn-1(i, k) represents the border of non-update area in (n-1)th 3-D view the i-th row kth layer
Residing columns, N represent the voxel columns of default neighborhood, herein N=3, Comn-1(i, j, k) represents the i-th row jth row kth layer
Base value is compensated, the fitting for compensating base value is fitted using linear mode.
5) the correction-compensation formula and compensation base value Com are utilizedn-1To the (n-1)th 3-D view Imgn-1Non- update area
Voxel value is modified compensation, and the voxel value obtained according to correction-compensation obtains the (n-1)th compensation image Fixn-1；The amendment is mended
Repaying formula is：
Wherein, Fixn-1(i, j, k) represents the voxel value of (n-1)th compensation image the i-th row jth row kth layer, JBn-1(i, k) table
Show the border of non-update area in the columns residing for (n-1)th 3-D view the i-th row kth layer, Imgn-1(i, j, k) represents the (n-1)th three
Tie up the voxel value of the i-th row jth row kth layer voxel in image, Comn-1(i, j, k) represents that (n-1)th 3-D view the i-th row jth arranges the
The compensation base value of k layers, Refm(i, j, k) represents the voxel value of m reference picture the i-th row jth row kth layers, and m represents current iteration
Number.
6) using parallel algebraic reconstruction SART technologies according to the (n-1)th compensation image Fixn-1To the data for projection of the n-th angle
ProjnRebuild, that is, obtain the n-th 3-D view Imgn；
To the (n-1)th compensation image Fixn-1Angularly n carries out forward projection, and bulb launches X ray under simulation current angular n
And penetrate on scanning project objects to detector, calculate each pixel of projected image (detector member) corresponding X ray penetrates road
Footpath, and will penetrate n-1 compensation image Fix is corresponded on path at this successivelyn-1Each voxel presses respective path weight on three dimensions
It is weighted cumulative, obtains the n-th estimated projection data Aprojn；
According to the n-th estimated projection data AProjnWith the data for projection Proj under the n-th anglenWeighted differences are carried out, will be estimated
The pixel value that projection corresponds to each pixel with acquired projections carries out difference processing, and as the X where the pixel (detector member)
The total length that ray penetrates scanning object is weighted processing to the pixel difference；
The result DProj that forward projection and weighted differences are obtainednBack projection's weighting amendment is carried out, it is angularly corresponding under n
Bulb, detector, the locus for scanning object, obtain projection difference result DProjnX ray corresponding to each coordinate, by the X
Ray path, which obtains, is pierced 3 d space coordinate corresponding to scanning object, by above-mentioned projection difference result by each on the path of place
The path weight value of individual voxel, weighted superposition to the (n-1)th compensation image Fixn-1Pixel value in corresponding three-dimensional coordinate pixel, is obtained
Update the n-th 3-D view Img of amendmentn。
7) judge whether n is equal to 21, if it is not, n then is updated into n+1, return to step 4) carry out successively under each angle not
Update area is calculated, 3-D view compensates, data for projection is rebuild, if so, the renewal for then completing all data for projection is rebuild, is entered
Step next step 8)；
8) by the n-th 3-D view ImgnAs m reconstruction image results Resm, judge whether m is equal to 5 (default iteration time
Number), if so, then using m reconstruction images result as final tomography composograph, if it is not, then to m reconstruction image results
ResmVoxel value and m reference image Rs efmVoxel value substitute into weighting more new formula, calculate m+1 reference image Rs efm+1；
Wherein, the weighting more new formula is：
Refm+1(i, j, k)=(Refm(i,j,k)+(m-1)pResm(i,j,k))/(1+(m-1)p),
Wherein, Refm+1(i, j, k) represents the voxel value of the i-th row jth row kth layer in m+1 reference pictures, Refm(i,j,
K) voxel value of the i-th row jth row kth layer in m reference pictures, Res are representedm(i, j, k) is represented the in m reconstruction image results
The voxel value of i row jth row kth layers, m represent current iteration number, and p is the renewal weight coefficient more than 0, and this application scenarios is set
For 1.
Fig. 6 show the DBT tomography composographs rebuild using the present invention program, and Fig. 7 is shown using traditional scheme
The DBT tomography composographs rebuild, as shown in fig. 7, using traditional SART method for reconstructing, because detector size is limited and
Angle is excessive, causes data for projection missing under some angles to be blocked, and causes the imaging FOV of the angle back projection missing area occur
Domain, in each projection renewal iterative operation, the mutation of absent region boundary voxel will be produced, more banding artifacts are formed, so as to influence
The effect that corresponding tomography composograph is rebuild.It will rebuild what is obtained using tomography composograph method for reconstructing provided by the invention
Tomography composograph (such as Fig. 6) is compared with the DBT tomography composographs shown in Fig. 7, it can be shown that utilizing the disconnected of the present invention
It is laminated into image rebuilding method, the reconstruction image of scanned object (measured target) and the structural information basic one of preferable die body
Cause；Occur the tomography composograph (such as Fig. 7) of more banding artifacts relative to border, pass through above-mentioned, tomography composograph of the invention
Tomography composograph obtained by method for reconstructing, it is clear that it is relatively sharp accurate, effectively solve above-mentioned wide-angle scanning and detector
Banding artifact problem caused by size-constrained, and peripheral components realization need not be increased, execution efficiency is high, and stability is strong.
With reference to figure 8, Fig. 8 show the tomography composograph reconstructing system structural representation of one embodiment, including：
Acquisition module 10, for obtaining the non-update area of the (n-1)th 3-D view, calculate m reference pictures and described the
The voxel difference of n-1 3-D views, interpolation fitting is carried out to non-update area according to voxel difference, obtains the non-update area
The compensation base value of each voxel；Wherein, n be more than 1 and less than or equal to n0 integer, n0 be transmitting X ray the child's hair twisted in a knot-childhood number of degrees, m
For current iteration number, initial value is set as 1；
Correcting module 20, for being carried out according to the compensation base value to the non-update area voxel value of the (n-1)th 3-D view
Correction-compensation, the voxel value obtained according to correction-compensation obtain the (n-1)th compensation image；
Module 30 is rebuild, for compensating projection number of the image to the n-th angle according to (n-1)th using parallel algebraic reconstruction technique
According to being rebuild, the n-th 3-D view is obtained；
Judge module 40, for judging whether n is equal to n0, if it is not, n then is updated into n+1, returns and perform acquisition module,
If so, then perform update module；
Update module 50, for using the n-th 3-D view as m reconstruction image results, judging whether m is equal to default change
Generation number, if it is not, setup module is performed, if so, then using m reconstruction images result as final tomography composograph；
Setup module 60, for m reconstruction images result and iterations m to be weighted into processing to m reference pictures,
M+1 reference pictures are obtained, and m is updated to m+1, returns and performs acquisition module.
In one embodiment, above-mentioned reconstruction module is further used for：
To the (n-1)th compensation image, angularly n carries out forward projection, obtains the n-th estimated projection data；
The weighted differences of pixel are carried out according to the data for projection under the n-th estimated projection data and the n-th angle；
After the result back projection processing that above-mentioned forward projection and weighted differences are obtained, to (n-1)th compensation each voxel of image
Amendment is weighted, the n-th 3-D view updated.
In one embodiment, above-mentioned setup module is further used for：
The voxel value of voxel value and m reference pictures to the m reconstruction image results substitutes into weighting more new formula, meter
Calculate the voxel value of m+1 reference pictures；Wherein, the weighting more new formula is：
Refm+1(i, j, k)=(Refm(i,j,k)+(m-1)pResm(i,j,k))/(1+(m-1)p),
Wherein, Refm+1(i, j, k) represents the voxel value of the i-th row jth row kth layer in m+1 reference pictures, Refm(i,j,
K) voxel value of the i-th row jth row kth layer in m reference pictures, Res are representedm(i, j, k) is represented the in m reconstruction image results
The voxel value of i row jth row kth layers, m represent current iteration number.
Tomography composograph reconstructing system provided by the invention and tomography composograph method for reconstructing one provided by the invention
One correspondence, the technical characteristic illustrated in the embodiment of the tomography composograph method for reconstructing and its advantage are applied to disconnected
In the laminated embodiment into image re-construction system, hereby give notice that.
Each technical characteristic of embodiment described above can be combined arbitrarily, to make description succinct, not to above-mentioned reality
Apply all possible combination of each technical characteristic in example to be all described, as long as however, the combination of these technical characteristics is not deposited
In contradiction, the scope that this specification is recorded all is considered to be.
Embodiment described above only expresses the several embodiments of the present invention, and its description is more specific and detailed, but simultaneously
Can not therefore it be construed as limiting the scope of the patent.It should be pointed out that come for one of ordinary skill in the art
Say, without departing from the inventive concept of the premise, various modifications and improvements can be made, these belong to the protection of the present invention
Scope.Therefore, the protection domain of patent of the present invention should be determined by the appended claims.
Claims (6)
1. a kind of tomography composograph method for reconstructing, it is characterised in that comprise the following steps：
S10, the non-update area for obtaining the (n-1)th 3-D view, calculate the body of m reference pictures and (n-1)th 3-D view
Plain difference, interpolation fitting is carried out to non-update area according to voxel difference, obtains the compensation base of non-each voxel of update area
Value；Wherein, n is the integer more than 1 and less than or equal to n0, and n0 is the child's hair twisted in a knot-childhood number of degrees of transmitting X ray, and m is current iteration number,
Initial value is set as 1；The data for projection of measured target is rebuild, obtains the initial reference image of measured target, i.e., first
Reference picture；In first time iterative process, the data for projection of first angle is rebuild, obtains the first 3-D view, profit
Relation is projected with light and scan geometry parameter determines the non-update area of the first 3-D view, according to first 3-D view
Non- update area, which calculates, accordingly obtains the first compensation image, and according to the data for projection under the compensation image and second angle
Carry out the reconstruction of the second 3-D view；
The calculating m reference pictures and the process of the voxel difference of (n-1)th 3-D view include：
The voxel value of the voxel value of the m reference pictures and the (n-1)th 3-D view is substituted into voxel differences formula respectively and calculates the
The voxel difference of n-1 3-D views；Wherein, the voxel differences formula is：
Diffn-1(i, j, k)=Imgn-1(i,j,k)-Refm(i,j,k)；
Wherein, Imgn-1(i, j, k) represents the voxel value of the i-th row jth row kth layer in the (n-1)th 3-D view, Refm(i, j, k) table
Show the voxel value of the i-th row jth row kth layer in m reference pictures, Diffn-1(i, j, k) represents the i-th row in the (n-1)th 3-D view
The voxel difference of jth row kth layer, m represent current iteration number；
It is described that interpolation fitting is carried out to non-update area according to voxel difference, obtain the compensation base of non-each voxel of update area
The process of value includes：
The voxel difference is substituted into the neighborhood correction value of the non-update area of neighboring mean value formula calculating；
The neighborhood correction value of the non-update area is substituted into the compensation base value of the non-update area of interpolation fitting formula calculating；
The neighboring mean value formula is：
<mrow>
<msub>
<mi>Ndiff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>N</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mrow>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
</mrow>
</munderover>
<msub>
<mi>Diff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>/</mo>
<mi>N</mi>
</mrow>
The interpolation fitting formula is：
<mrow>
<msub>
<mi>Com</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>Ndiff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mfrac>
<mrow>
<mi>j</mi>
<mo>-</mo>
<mi>J</mi>
</mrow>
<mrow>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>J</mi>
</mrow>
</mfrac>
<mo>,</mo>
<mi>j</mi>
<mo>&Element;</mo>
<mo>(</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
<mo>)</mo>
<mo>,</mo>
<mi>J</mi>
<mo>&rsqb;</mo>
<mo>,</mo>
</mrow>
Wherein, Ndiffn-1(i, k) represents the neighborhood correction value of the i-th row kth layer in the (n-1)th 3-D view, and J represents that (n-1)th is three-dimensional
The border columns of voxel, JB in imagen-1(i, k) represents the border of non-update area in (n-1)th 3-D view the i-th row kth layer institute
The columns at place, N represent the voxel columns of default neighborhood, Comn-1(i, j, k) represents (n-1)th 3-D view the i-th row jth row kth layer
Compensation base value；
S20, compensation is modified to the non-update area voxel value of the (n-1)th 3-D view according to the compensation base value, according to repairing
Just compensate obtained voxel value and obtain the (n-1)th compensation image；
S30, using parallel algebraic reconstruction technique according to (n-1)th compensation image the data for projection of the n-th angle is rebuild, obtain
N-th 3-D view；
S40, judge whether n is equal to n0, if it is not, n then is updated into n+1, return to step S10, if so, then entering step S50；
S50, using the n-th 3-D view as m reconstruction image results, judge whether m is equal to default iterations, if it is not, entering
Enter step S60, if so, then using m reconstruction images result as final tomography composograph；
S60, processing is weighted to m reference pictures according to m reconstruction images result and iterations m, obtains m+1 references
Image, and m is updated to m+1, return to step S10；The step S60 includes：
The voxel value of voxel value and m reference pictures to the m reconstruction image results substitutes into weighting more new formula, calculates the
The voxel value of m+1 reference pictures；Wherein, the weighting more new formula is：
Refm+1(i, j, k)=(Refm(i,j,k)+(m-1)pResm(i,j,k))/(1+(m-1)p),
Wherein, Refm+1(i, j, k) represents the voxel value of the i-th row jth row kth layer in m+1 reference pictures, Refm(i, j, k) table
Show the voxel value of the i-th row jth row kth layer in m reference pictures, Resm(i, j, k) represents the i-th row in m reconstruction image results
The voxel value of jth row kth layer, m represent current iteration number, and p is renewal weight coefficient.
2. tomography composograph method for reconstructing according to claim 1, it is characterised in that also wrapped before the step S10
Include：
Backprojection reconstruction is carried out to multiple data for projection that detector receives, the 3-D view that backprojection reconstruction is obtained is arranged to
First reference picture；
The data for projection of first angle is rebuild using parallel algebraic reconstruction technique, obtains the first 3-D view；
The non-update area of the first 3-D view is obtained, calculates the voxel differences of the first reference picture and first 3-D view
Value, interpolation fitting is carried out to non-update area according to voxel difference, obtains the compensation base value of non-each voxel of update area；
The compensation base value and the first reference picture are modified compensation to the voxel value of the first 3-D view, mended according to amendment
The voxel value for repaying to obtain obtains the first compensation image；
The data for projection of second angle is rebuild according to the first compensation image using parallel algebraic reconstruction technique, obtains second
3-D view.
3. tomography composograph method for reconstructing according to claim 1, it is characterised in that the step S30 includes：
To the (n-1)th compensation image, angularly n carries out forward projection, obtains the n-th estimated projection data；
The weighted differences of pixel are carried out according to the data for projection under the n-th estimated projection data and the n-th angle；
After the result back projection processing that above-mentioned forward projection and weighted differences are obtained, (n-1)th each voxel of compensation image is carried out
Weighting amendment, the n-th 3-D view updated.
4. tomography composograph method for reconstructing according to claim 1, it is characterised in that it is described according to compensation base value to the
The process that the voxel value of n-1 3-D views is modified compensation includes：
The voxel value of the compensation base value and the (n-1)th 3-D view is substituted into correction-compensation formula respectively and carries out computing；It is described to repair
Positive compensation formula is：
<mrow>
<msub>
<mi>Fix</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<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>
<msub>
<mi>Img</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>j</mi>
<mo>&Element;</mo>
<mo>&lsqb;</mo>
<mn>0</mn>
<mo>,</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>Ref</mi>
<mi>m</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msub>
<mi>Com</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>j</mi>
<mo>&Element;</mo>
<mrow>
<mo>(</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>J</mi>
<mo>&rsqb;</mo>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>,</mo>
</mrow>
Wherein, Fixn-1(i, j, k) represents the voxel value of (n-1)th compensation image the i-th row jth row kth layer, JBn-1(i, k) is represented not
The border of update area is in the columns residing for (n-1)th 3-D view the i-th row kth layer, Imgn-1(i, j, k) represents the (n-1)th graphics
The voxel value of the i-th row jth row kth layer voxel, Com as inn-1(i, j, k) represents (n-1)th 3-D view the i-th row jth row kth layer
Compensation base value, Refm(i, j, k) represents the voxel value of m reference picture the i-th row jth row kth layers.
A kind of 5. tomography composograph reconstructing system, it is characterised in that including：
Acquisition module, for obtaining the non-update area of the (n-1)th 3-D view, calculate m reference pictures and the described (n-1)th three-dimensional
The voxel difference of image, interpolation fitting is carried out to non-update area according to voxel difference, obtains non-each voxel of update area
Compensation base value；Wherein, n is the integer more than 1 and less than or equal to n0, and n0 is the child's hair twisted in a knot-childhood number of degrees of transmitting X ray, and m is current
Iterations, initial value are set as 1；The data for projection of measured target is rebuild, obtains the initial reference figure of measured target
Picture, i.e. the first reference picture；In first time iterative process, the data for projection of first angle is rebuild, obtains the one or three
Image is tieed up, the non-update area that relation and scan geometry parameter determine the first 3-D view is projected using light, according to described the
One 3-D view and corresponding non-update area, which calculate, obtains the first compensation image, and according under the compensation image and second angle
Data for projection carry out the second 3-D view reconstruction；
The calculating m reference pictures and the process of the voxel difference of (n-1)th 3-D view include：
The voxel value of the voxel value of the m reference pictures and the (n-1)th 3-D view is substituted into voxel differences formula respectively and calculates the
The voxel difference of n-1 3-D views；Wherein, the voxel differences formula is：
Diffn-1(i, j, k)=Imgn-1(i,j,k)-Refm(i,j,k)；
Wherein, Imgn-1(i, j, k) represents the voxel value of the i-th row jth row kth layer in the (n-1)th 3-D view, Refm(i, j, k) table
Show the voxel value of the i-th row jth row kth layer in m reference pictures, Diffn-1(i, j, k) represents the i-th row in the (n-1)th 3-D view
The voxel difference of jth row kth layer, m represent current iteration number；
It is described that interpolation fitting is carried out to non-update area according to voxel difference, obtain the compensation base of non-each voxel of update area
The process of value includes：
The voxel difference is substituted into the neighborhood correction value of the non-update area of neighboring mean value formula calculating；
The neighborhood correction value of the non-update area is substituted into the compensation base value of the non-update area of interpolation fitting formula calculating；
The neighboring mean value formula is：
<mrow>
<msub>
<mi>Ndiff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>j</mi>
<mo>=</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>N</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
<mrow>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
</mrow>
</munderover>
<msub>
<mi>Diff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>/</mo>
<mi>N</mi>
</mrow>
The interpolation fitting formula is：
<mrow>
<msub>
<mi>Com</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>j</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>Ndiff</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mfrac>
<mrow>
<mi>j</mi>
<mo>-</mo>
<mi>J</mi>
</mrow>
<mrow>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mi>J</mi>
</mrow>
</mfrac>
<mo>,</mo>
<mi>j</mi>
<mo>&Element;</mo>
<mo>(</mo>
<msub>
<mi>JB</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>(</mo>
<mrow>
<mi>i</mi>
<mo>,</mo>
<mi>k</mi>
</mrow>
<mo>)</mo>
<mo>,</mo>
<mi>J</mi>
<mo>&rsqb;</mo>
<mo>,</mo>
</mrow>
Wherein, Ndiffn-1(i, k) represents the neighborhood correction value of the i-th row kth layer in the (n-1)th 3-D view, and J represents that (n-1)th is three-dimensional
The border columns of voxel, JB in imagen-1(i, k) represents the border of non-update area in (n-1)th 3-D view the i-th row kth layer institute
The columns at place, N represent the voxel columns of default neighborhood, Comn-1(i, j, k) represents (n-1)th 3-D view the i-th row jth row kth layer
Compensation base value；
Correcting module, for being modified benefit to the non-update area voxel value of the (n-1)th 3-D view according to the compensation base value
Repay, the voxel value obtained according to correction-compensation obtains the (n-1)th compensation image；
Module is rebuild, for being carried out using parallel algebraic reconstruction technique according to the (n-1)th compensation image to the data for projection of the n-th angle
Rebuild, obtain the n-th 3-D view；
Judge module, for judging whether n is equal to n0, if it is not, n then is updated into n+1, returns and perform acquisition module, if so, then
Perform update module；
Update module, for using the n-th 3-D view as m reconstruction image results, judging whether m is equal to default iteration time
Number, if it is not, setup module is performed, if so, then using m reconstruction images result as final tomography composograph；
Setup module, for m reconstruction images result and iterations m to be weighted into processing to m reference pictures, obtain
M+1 reference pictures, and m is updated to m+1, return and perform acquisition module；The setup module is further used for：
The voxel value of voxel value and m reference pictures to the m reconstruction image results substitutes into weighting more new formula, calculates the
The voxel value of m+1 reference pictures；Wherein, the weighting more new formula is：
Refm+1(i, j, k)=(Refm(i,j,k)+(m-1)pResm(i,j,k))/1+(m-1)p,
Wherein, Refm+1(i, j, k) represents the voxel value of the i-th row jth row kth layer in m+1 reference pictures, Refm(i, j, k) table
Show the voxel value of the i-th row jth row kth layer in m reference pictures, Resm(i, j, k) represents the i-th row in m reconstruction image results
The voxel value of jth row kth layer, m represent current iteration number.
6. tomography composograph reconstructing system according to claim 5, it is characterised in that the reconstruction module is further used
In：
To the (n-1)th compensation image, angularly n carries out forward projection, obtains the n-th estimated projection data；
The weighted differences of pixel are carried out according to the data for projection under the n-th estimated projection data and the n-th angle；
After the result back projection processing that above-mentioned forward projection and weighted differences are obtained, (n-1)th each voxel of compensation image is carried out
Weighting amendment, the n-th 3-D view updated.
Priority Applications (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610796748.6A CN106373165B (en) | 2016-08-31 | 2016-08-31 | Tomography composograph method for reconstructing and system |
| PCT/CN2016/098738 WO2018040126A1 (en) | 2016-08-31 | 2016-09-12 | Tomosynthesis image reconstruction method and system |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610796748.6A CN106373165B (en) | 2016-08-31 | 2016-08-31 | Tomography composograph method for reconstructing and system |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN106373165A CN106373165A (en) | 2017-02-01 |
| CN106373165B true CN106373165B (en) | 2017-11-28 |
Family
ID=57899881
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201610796748.6A Active CN106373165B (en) | 2016-08-31 | 2016-08-31 | Tomography composograph method for reconstructing and system |
Country Status (2)
| Country | Link |
|---|---|
| CN (1) | CN106373165B (en) |
| WO (1) | WO2018040126A1 (en) |
Families Citing this family (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108010096A (en) | 2017-11-30 | 2018-05-08 | 上海联影医疗科技有限公司 | CBCT image rebuilding methods, device and CBCT equipment |
| DE102017223603B4 (en) * | 2017-12-21 | 2020-01-09 | Siemens Healthcare Gmbh | Process for the reconstruction of a three-dimensional image data set, recorded with a biplan x-ray device, biplan x-ray device, computer program and electronically readable data carrier |
| US20210233293A1 (en) * | 2018-05-04 | 2021-07-29 | Our United Corporation | Low-dose imaging method and apparatus |
| CN109147035B (en) * | 2018-08-03 | 2023-08-08 | 上海电气集团股份有限公司 | Three-dimensional model display method and system |
| EP3660790A1 (en) * | 2018-11-28 | 2020-06-03 | Koninklijke Philips N.V. | System for reconstructing an image of an object |
| CN109658465B (en) * | 2018-12-07 | 2023-07-04 | 广州华端科技有限公司 | Data processing in image reconstruction process, image reconstruction method and device |
| CN110335204B (en) * | 2019-05-07 | 2021-06-01 | 中国人民解放军陆军工程大学 | Thermal imaging image enhancement method |
| CN110706338B (en) * | 2019-09-30 | 2023-05-02 | 东软医疗系统股份有限公司 | Image reconstruction method, device, CT equipment and CT system |
| CN112258627B (en) * | 2020-09-18 | 2023-09-15 | 中国科学院计算技术研究所 | Local fault three-dimensional reconstruction system |
| CN115147378B (en) * | 2022-07-05 | 2023-07-25 | 哈尔滨医科大学 | CT image analysis and extraction method |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7672424B2 (en) * | 2004-10-08 | 2010-03-02 | Koninklijke Philips Electronics N.V. | Image reconstruction with voxel dependent interpolation |
| CN102419866B (en) * | 2010-09-27 | 2013-08-21 | 北京农业智能装备技术研究中心 | Construction method of filter function for CT (computerized tomography) image tomography reconstruction and tomography reconstruction method |
| CN102103757B (en) * | 2010-12-27 | 2012-09-19 | 中国科学院深圳先进技术研究院 | Cone beam image rebuilding method and device |
| CN103971349B (en) * | 2013-01-30 | 2017-08-08 | 上海西门子医疗器械有限公司 | Computed tomography images method for reconstructing and ct apparatus |
| CN103455989B (en) * | 2013-09-24 | 2016-06-15 | 南京大学 | A kind of method improving limited angle CT image quality in conjunction with ultrasonoscopy |
| US9600924B2 (en) * | 2014-02-05 | 2017-03-21 | Siemens Aktiengesellschaft | Iterative reconstruction of image data in CT |
| CN105488823B (en) * | 2014-09-16 | 2019-10-18 | 株式会社日立制作所 | CT image rebuilding method, CT equipment for reconstructing image and CT system |
| CN105184835B (en) * | 2015-09-15 | 2018-11-06 | 上海联影医疗科技有限公司 | Mammary gland tomographic image reconstructing process and device |
-
2016
- 2016-08-31 CN CN201610796748.6A patent/CN106373165B/en active Active
- 2016-09-12 WO PCT/CN2016/098738 patent/WO2018040126A1/en active Application Filing
Also Published As
| Publication number | Publication date |
|---|---|
| WO2018040126A1 (en) | 2018-03-08 |
| CN106373165A (en) | 2017-02-01 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN106373165B (en) | Tomography composograph method for reconstructing and system | |
| KR101728046B1 (en) | Tomography apparatus and method for reconstructing a tomography image thereof | |
| US5552605A (en) | Motion correction based on reprojection data | |
| RU2510080C2 (en) | Image processing device, image processing method and long-term information storage medium | |
| US20100172472A1 (en) | Collecting images for image stitching with rotating a radiation detector | |
| US9848844B2 (en) | Iterative reconstruction process | |
| JP2011503570A (en) | Apparatus and method for forming an attenuation map | |
| US11179128B2 (en) | Methods and systems for motion detection in positron emission tomography | |
| US20140212018A1 (en) | System optics in at least in one of backprojection and forward projection for model-based iterative reconstruction | |
| CN104274201A (en) | Method, system and equipment for tomography of mammary gland and image acquisition and processing method | |
| CN103797517A (en) | Method of image reconstruction for a filtered back projection in limited angle tomography | |
| US9965875B2 (en) | Virtual projection image method | |
| CN109146800B (en) | Cone beam computed tomography method for correcting image and system | |
| CN108992083B (en) | Cone beam computed tomography method for correcting image and system | |
| US10062168B2 (en) | 5D cone beam CT using deformable registration | |
| CN102062740B (en) | Cone-beam CT (Computed Tomography) scanning imaging method and system | |
| US10089758B1 (en) | Volume image reconstruction using projection decomposition | |
| US9237873B2 (en) | Methods and systems for CT projection domain extrapolation | |
| CN111956248A (en) | X-ray imaging method, device, equipment and storage medium | |
| US10966670B2 (en) | Imaging system and method for dual-energy and computed tomography | |
| CN105631908A (en) | PET image reconstruction method and device | |
| CN108024779A (en) | Computed tomography image generating means | |
| CN109620273A (en) | A kind of quick CBCT algorithm for reconstructing calculating short scanning weight in real time | |
| CN101810491B (en) | Reduction of artifacts caused by movement of an X-ray tube in object reconstruction | |
| CN108366772B (en) | Image processing apparatus and method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | 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 |