Electrical property tomography method, system, device and medium based on instantaneous linearization
By measuring the magnetic resonance radio frequency field, calculating the point data information of the transmission field, combining Maxwell's equations and Newton's iteration method, and adding constraints, the error problem caused by the assumption of uniform distribution of electrical characteristics in the existing technology is solved, and higher-precision electrical characteristic tomography is achieved.
Patent Information
- Application Number
- CN202211343568.4
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2022-10-31
- Publication Date
- 2025-09-26
- Estimated Expiration
- 2042-10-31
AI Technical Summary
In existing electrical property tomography methods, the uniformity assumption leads to large calculation errors, and the Helmholtz equation is difficult to accurately reconstruct the non-uniform electrical property distribution, affecting imaging accuracy.
An electrical property tomography method based on instantaneous linearization is adopted. By measuring the magnetic resonance radio frequency field, the point data information of the transmitting field is calculated. Combining Maxwell's equations and Newton's iteration method, the constraints K, R1 and R2 are added to construct an iterative equation, and the permittivity and conductivity are calculated step by step.
The precision and accuracy of electrical property tomography are improved, the model simplification is reduced, and the electrical property distribution is more consistent with the actual situation.
Smart Images

Figure CN115808650B_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the field of nuclear magnetic resonance imaging technology, and in particular relates to an electrical property tomography method, system, equipment and medium based on instantaneous linearization. Background Art
[0002] In high-field MRI systems, the frequency increases with increasing field strength, while the wavelength decreases until it matches the wavelength of the object being scanned. This state induces resonance, strengthening the interaction between the object being measured and the RF field, a phenomenon that can be exploited to achieve imaging. This high-field MRI-based method, called MR-ERT, is a new quantitative MRI imaging method.
[0003] Biological tissues contain large quantities of water and small amounts of inorganic salts, substances with strong electrical properties. The distribution of these substances varies between tissues. This characteristic makes biological tissue a highly electrically charged object, which can be exploited for imaging. Unlike traditional tomographic imaging methods, dielectric tomography is a quantitative measurement, allowing for quantitative analysis. This offers inherent advantages for addressing problems such as early cancer analysis, which are difficult to address with traditional tomographic imaging. Furthermore, dielectric properties influence the propagation of electrical signals within the brain and are highly correlated with information such as brain tissue type, activity, function, aging, and pathological changes. Therefore, MR-EPT has the potential to provide a novel research tool for studying brain cognition and brain disease. The distribution of dielectric properties in the human body also forms the basis for calculating the specific absorption rate (SAR), which can be used to assess human safety in microwave and radiofrequency fields. Thermal effects have always been a serious safety concern for ultra-high-field MRI systems. Accurate imaging of the dielectric properties of human tissue can facilitate SAR research and promote the safer application of ultra-high-field MRI systems.
[0004] Dielectric tomography involves two steps. The first step is to acquire radiofrequency field information. The second step is to derive the equations for dielectric tomography based on Maxwell's equations. A mathematical model is then constructed to solve the electrical property distribution. The model uses the acquired radiofrequency field information as input to reconstruct the electrical property distribution. This second step is a key research topic for electrical tomography algorithms: how to construct a correct mathematical model to more accurately reconstruct the electrical property distribution.
[0005] Currently, most mainstream algorithms for electrical property reconstruction ignore minor issues or simplify the complex parameters of the equations themselves. One classic algorithm assumes that the electrical properties themselves are uniformly distributed. This method is called the Helmholtz EPT. The drawback of this method is very obvious: the actual distribution of electrical properties cannot be completely uniform, which means that the premise of the method itself does not meet realistic conditions. This will invalidate the uniformity assumption, leading to serious calculation errors and generating a lot of erroneous data in non-uniform areas. Moreover, the Helmholtz equation itself is a high-order, nonlinear, non-homogeneous partial differential complex domain equation, which means that without simplification, the equation itself is difficult to solve. Summary of the Invention
[0006] The object of the present invention is to provide an electrical property tomography method based on instantaneous linearization, which can obtain electrical property tomography more accurately while minimizing the simplification or modification of the problem model.
[0007] The present invention is achieved through the following technical solutions:
[0008] An electrical property tomography method based on instantaneous linearization comprises the following steps:
[0009] S1. Measure the magnetic resonance radio frequency field and calculate the transmit field based on the radio frequency field;
[0010] S2. Calculate the point data information of the launch field using a numerical method. The point data information includes the first-order derivative and the second-order derivative along the x-axis, y-axis, and z-axis at each point in space of the launch field.
[0011] S3. Based on the launch site and the point data information of the launch site, according to Maxwell's equations, calculate the initial value of the permittivity and the initial value of the conductivity at each point in space, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space;
[0012] S4. Based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, the iterative equation is obtained according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1):
[0013]
[0014] Where u n represents u of the nth generation, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector;
[0015] S5. Based on the initial values of the permittivity and conductivity at each point in the space, determine the initial values of the iterative equation, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations;
[0016] S6. Output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
[0017] Furthermore, in the iterative equation, f(u) is as shown in formula (2):
[0018] f(u)=Au+Cu 2 (2)
[0019] Where u 2 is a column vector, u 2 The general formula for each term is A is a matrix, as shown in formula (3):
[0020]
[0021] In the formula, the variable As shown in formula (4), the variable γ is shown in formula (5):
[0022]
[0023]
[0024] C is a matrix, as shown in formula (6):
[0025]
[0026] Furthermore, in the iterative equation, the damping factor vector α is calculated as follows:
[0027] Build an n x * y * z The matrix Q of , the value of each point of Q is shown in formula (7):
[0028]
[0029] Where σ m and ε m All are parameters;
[0030] Suppose there is a matrix V. For each point v in V, there is a corresponding m*m matrix T in Q with the coordinate of v as the center. The value of v is shown in formula (23):
[0031] v=min(T) (8)
[0032] The damping factor vector α is obtained by mapping the matrix V one by one to the damping factor vector α.
[0033] Furthermore, the steps of determining the constraints K, R1, and R2 from the initial values of the permittivity and the conductivity at each point in the space include:
[0034] For the initial value of the permittivity and the initial value of the conductivity at each point in the space, determine whether the initial value of the permittivity and the initial value of the conductivity along all directions are 0;
[0035] If so, the initial values of permittivity and conductivity are recorded as the target initial values of permittivity and conductivity, and the target initial values of permittivity and conductivity do not participate in the iteration. The constraint term K is constructed based on the target initial values of permittivity and conductivity.
[0036] If not, the initial values of permittivity and conductivity are recorded as the initial values of iterative permittivity and iterative conductivity;
[0037] The constraint R1 is of size 1*(n x *n y *n z ), the general formula of each item in the constraint term R1 is shown in formula (9):
[0038]
[0039] Where λ R1 is the coefficient of constraint R1, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, and 0 and (x0, y0, z0) satisfy the relationship shown in formula (10):
[0040] n0=n z *n y *x0+n z *y0+z0 (10)
[0041] f R1 As shown in formula (11):
[0042]
[0043] Where θ r is a parameter;
[0044] Constant term θ t0 As shown in formula (12):
[0045]
[0046] Where, the constant p is the number of target permittivity initial values;
[0047] The constraint R2 is a scale of 1*(nx *n y *n z ), the constraint term R2 is shown in formula (13):
[0048]
[0049] Where λ R2 is the coefficient of constraint R2, is the Laplace operator, u is a scale of 1*(n x *n y *n z ), each term of u is expressed as follows:
[0050]
[0051] Furthermore, based on the initial values of the permittivity and the conductivity at each point in the space, the initial values of the iterative equation are determined, the initial values are substituted into the iterative equation for iterative calculation, and the permittivity and conductivity at each point in the space are obtained after several iterations. The steps include:
[0052] S51. Preprocess the initial value of the iterative permittivity and the initial value of the iterative conductivity. Based on the preprocessed initial value of the iterative permittivity and the initial value of the iterative conductivity, calculate the initial value u0 of the iterative equation according to formula (15):
[0053] u0=σ0+iωε0 (15)
[0054] Where, σ0 is the initial value of iterative conductivity, ε0 is the initial value of iterative permittivity;
[0055] S52, substituting the initial value u0 into the iterative equation to obtain a new initial value, and determining whether the average value of each item in the new initial value is less than a preset iteration threshold;
[0056] S53, if not, repeat step S52;
[0057] S54. If yes, the permittivity and conductivity at each point in the space are obtained according to the new initial value after iteration, the target permittivity initial value, and the target conductivity initial value.
[0058] Furthermore, the step of calculating the point data information of the launch site using a numerical method includes:
[0059] According to formula (16), the first-order derivative of the launch field along the x-axis at each point in space is calculated as follows:
[0060]
[0061] According to formula (17), the first-order derivative of the launch field along the y-axis at each point in space is calculated as follows:
[0062]
[0063] According to formula (18), the first-order derivative of the launch field along the z-axis at each point in space is calculated as follows:
[0064]
[0065] According to formula (19), the second-order derivative of the launch field along the x-axis at each point in space is calculated as
[0066]
[0067] According to formula (20), the second-order derivative of the launch field along the y-axis at each point in space is calculated as
[0068]
[0069] According to formula (21), the second-order derivative of the launch field along the z-axis at each point in space is calculated as
[0070]
[0071] Among them, dx, dy, and dz are the lengths of the unit grid in the x-axis, y-axis, and z-axis directions after spatial discretization, respectively.
[0072] Furthermore, the steps of measuring the magnetic resonance radio frequency field and calculating the transmit field based on the radio frequency field include:
[0073] Measure the components B of the magnetic resonance radio frequency field along the x-axis, y-axis, and z-axis of the spatial coordinate system x 、B y 、B z ;
[0074] Calculate the launch field B according to formula (22):
[0075]
[0076] Where i is a complex unit.
[0077] The present invention also provides an electrical property tomography system based on instantaneous linearization, comprising:
[0078] A measurement module is used to measure the magnetic resonance radio frequency field and calculate the transmission field based on the radio frequency field;
[0079] A first calculation module is used to calculate point data information of the launch field using a numerical method, where the point data information includes first-order derivatives and second-order derivatives of the launch field along the x-axis, y-axis, and z-axis at each point in space;
[0080] The second calculation module is used to calculate the initial value of the permittivity and the initial value of the conductivity at each point in space based on the launch site and the point data information of the launch site according to Maxwell's equations, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space;
[0081] The module is used to obtain the iterative equation based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1):
[0082]
[0083] Where u n represents u of the nth generation, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector;
[0084] An iterative module is used to determine the initial values of the iterative equation based on the initial values of the permittivity and conductivity at each point in the space, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations;
[0085] The output module is used to output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
[0086] The present invention also discloses an electronic device, comprising a memory and a processor, wherein the memory stores a computer program, and the processor implements the steps of any of the above methods when executing the computer program.
[0087] The present invention also discloses a computer-readable storage medium on which a computer program is stored. When the computer program is executed by a processor, the steps of any of the above methods are implemented.
[0088] Compared with the existing technology, the beneficial effects of the present invention are: the high-order nonlinear nonhomogeneous partial differential equations derived from the original Maxwell equations and the Helmholtz equations are instantaneously linearized through the improved Newton iteration method to obtain iterative equations, and constraint terms are added according to the physical model and actual conditions to form a complete iterative equation with initial values, thereby obtaining an electrical characteristic distribution solution with higher accuracy and more in line with the actual situation, thereby being able to obtain electrical characteristic tomography more accurately while minimizing the simplification or change of the problem model. BRIEF DESCRIPTION OF THE DRAWINGS
[0089] Figure 1Flow chart of the steps of the instantaneous linear electrical property tomography method of the present invention;
[0090] Figure 2 Comparison of the reconstructed conductivity (top) and permittivity (bottom) images of the Duke human brain tumor simulation model using the instantaneous linearization-based electrical property tomography method of the present invention, the Helmholtz EPT method, and the modified crEPT method.
[0091] Figure 3 The results of reconstructing conductivity (top) and permittivity (bottom) images using the improved crEPT method based on the instantaneous linearized electrical property tomography method and the Helmholtz EPT method on another Duke human brain tumor simulation model are shown in the figure.
[0092] Figure 4 Schematic diagram of the modules of the instantaneous linearization-based electrical property tomography system of the present invention;
[0093] Figure 5 This is a schematic structural diagram of an electronic device according to an embodiment of the present invention;
[0094] Figure 6 This is a schematic block diagram of the structure of an embodiment of a computer-readable storage medium of the present invention. DETAILED DESCRIPTION
[0095] To make the objectives, technical solutions, and advantages of the embodiments of the present invention more clear, the technical solutions of the embodiments of the present invention will be clearly and completely described below in conjunction with the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Generally, the components of the embodiments of the present invention described and shown in the drawings herein can be arranged and designed in various different configurations.
[0096] Therefore, the following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention as claimed, but rather merely represents selected embodiments of the present invention. All other embodiments derived by persons of ordinary skill in the art based on the embodiments of the present invention without creative effort shall fall within the scope of protection of the present invention.
[0097] It should be noted that similar reference numerals and letters represent similar items in the following drawings. Therefore, once an item is defined in one drawing, it does not need to be further defined or explained in subsequent drawings. At the same time, in the description of the present invention, the terms "first", "second", etc. are used only to distinguish the description and should not be understood as indicating or implying relative importance.
[0098] It should be noted that, in this document, relational terms such as first and second, etc., are used only to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply the existence of any such actual relationship or order between these entities or operations. Moreover, the terms "comprises," "comprising," or any other variants thereof are intended to cover non-exclusive inclusion, so that a process, method, article, or device comprising a series of elements includes not only those elements, but also other elements not explicitly listed, or elements inherent to such process, method, article, or device. In the absence of further limitations, an element defined by the phrase "comprising a ..." does not exclude the presence of other identical elements in the process, method, article, or device comprising the element.
[0099] In the description of the present invention, it should be noted that the terms "upper", "lower", "inside", "outside", etc. indicate orientations or positional relationships based on the orientations or positional relationships shown in the accompanying drawings, or are the orientations or positional relationships in which the inventive product is usually placed when in use. They are only for the convenience of describing the present invention and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation. Therefore, they should not be understood as limiting the present invention.
[0100] See also Figure 1 , Figure 1 The present invention is a flowchart of the steps of the instantaneous linearization electrical property tomography method. The instantaneous linearization electrical property tomography method includes the following steps:
[0101] S1. Measure the magnetic resonance radio frequency field and calculate the transmit field based on the radio frequency field;
[0102] S2. Calculate the point data information of the launch field using a numerical method. The point data information includes the first-order derivative and the second-order derivative along the x-axis, y-axis, and z-axis at each point in space of the launch field.
[0103] S3. Based on the launch site and the point data information of the launch site, according to Maxwell's equations, calculate the initial value of the permittivity and the initial value of the conductivity at each point in space, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space;
[0104] S4. Based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, the iterative equation is obtained according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1):
[0105]
[0106] Where u nrepresents u of the nth generation, i.e., the iterative value obtained after the nth iteration, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector;
[0107] S5. Based on the initial values of the permittivity and conductivity at each point in the space, determine the initial values of the iterative equation, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations;
[0108] S6. Output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
[0109] In step S1 above, the MRI device, in conjunction with the orthogonal coil and a specific B1 scanning sequence, can obtain the radio frequency field B in the orthogonal coil measurement space, that is, the magnetic resonance radio frequency field B. Then, based on the obtained radio frequency field B, the transmit field B is calculated.
[0110] Furthermore, in step S1, the steps of measuring the magnetic resonance radio frequency field and calculating the transmit field based on the radio frequency field include:
[0111] S11, measure the components B of the magnetic resonance radio frequency field along the x-axis, y-axis and z-axis directions of the spatial coordinate system x 、B y 、B z ;
[0112] S12. Calculate the launch site B according to formula (22):
[0113]
[0114] Where i is a complex unit.
[0115] In the above steps S11 to S12, the unique vector of the radio frequency field B measured by the instrument is (B x ,B y ,B z ), where B x 、B y 、B z are the components of the radio frequency field B along the x-axis, y-axis and z-axis of the spatial coordinate system. The transmitting field B is given by the formula calculate.
[0116] In the above step S2, the numerical method can adopt the finite difference method. Based on the launch field B, the first-order derivative and the second-order derivative of the launch field along the x-axis, y-axis and z-axis at each point in space are calculated by the finite difference method.
[0117] Furthermore, the step of calculating the point data information of the launch site using a numerical method includes:
[0118] According to formula (16), the first-order derivative of the launch field along the x-axis at each point in space is calculated as follows:
[0119]
[0120] According to formula (17), the first-order derivative of the launch field along the y-axis at each point in space is calculated as follows:
[0121]
[0122] According to formula (18), the first-order derivative of the launch field along the z-axis at each point in space is calculated as follows:
[0123]
[0124] According to formula (19), the second-order derivative of the launch field along the x-axis at each point in space is calculated as
[0125]
[0126] According to formula (20), the second-order derivative of the launch field along the y-axis at each point in space is calculated as
[0127]
[0128] According to formula (21), the second-order derivative of the launch field along the z-axis at each point in space is calculated as
[0129]
[0130] Where dx, dy, and dz are the lengths of the unit grid in the x-axis, y-axis, and z-axis directions after spatial discretization, respectively, and i is a subscript indicating a specific point in space.
[0131] In the above step S3, based on the launch field and the point data information of the launch field, according to the principle of MREPT (magnetic resonance electrical properties, magnetic resonance tomography of electrical properties of human tissue), the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at each point in the space are calculated, and then the calculated initial value of the permittivity ε0 and the initial value of the conductivity σ0 at each point in the space are processed according to the preset processing rules to determine the constraints K, R1 and R2.
[0132] Furthermore, in step S3, the step of determining the constraints K, R1, and R2 from the initial values of the permittivity and the conductivity at each point in the space includes:
[0133] S31. For the initial value of the permittivity and the initial value of the conductivity at each point in the space, determine whether the initial values of the permittivity and the initial values of the conductivity along all directions are 0;
[0134] S32. If yes, the initial permittivity value and the initial conductivity value are recorded as the target permittivity value and the target conductivity value, and the target permittivity value and the target conductivity value are not involved in the iteration. The constraint term K is constructed based on the target permittivity value and the target conductivity value.
[0135] S33. If not, the initial value of the permittivity and the initial value of the conductivity are recorded as the iterative initial value of the permittivity and the iterative initial value of the conductivity;
[0136] S34, the constraint R1 is a scale of 1*(n x *n y *n z ), the general formula of each item in the constraint term R1 is shown in formula (9):
[0137]
[0138] Where λ R1 is the coefficient of constraint R1, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, and n0 and (x0, y0, z0) satisfy the relationship shown in formula (10):
[0139] n0=n z *n y *x0+n z *y0+z0 (10)
[0140] f R1 As shown in formula (11):
[0141]
[0142] Where θ r It is a parameter, and its value is usually 0.2;
[0143] Constant term θ t0 As shown in formula (12):
[0144]
[0145] Where, the constant p is the number of target permittivity initial values;
[0146] S35, the constraint R2 is a scale of 1*(n x *n y *n z ), the constraint term R2 is shown in formula (13):
[0147]
[0148] Where λ R2 is the coefficient of constraint R2, is the Laplace operator, u is a scale of 1*(n x *n y *n z ), each term of u is expressed as follows:
[0149]
[0150] In steps S31 to S33 above, referring to the Helmholtz EPT method, it is assumed that the electrical characteristic parameters are uniformly distributed in space, that is, the derivatives of the electrical characteristic parameters along all directions are 0 at every point in space. This allows the partial differential equation to be converted into a linear equation, facilitating the acquisition of initial values for iteration, where the electrical characteristic parameters include permittivity and conductivity. In the real spatial distribution, some initial permittivity values ε0 and conductivity values σ0 do satisfy the assumption that the derivatives along all directions are equal to 0. These initial permittivity values ε0 and conductivity values σ0 are already very accurate and do not need to be solved again. Instead, the iterative equation itself can be constrained by converting the unsolved quantities into known quantities. This reduces the computational complexity of subsequent iterative calculations, thereby reducing the amount of computation required. Therefore, based on the assumptions of the Helmholtz EPT method, for the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at each point in space, analyze whether the initial value of the permittivity ε0 and the initial value of the conductivity σ0 along all directions are 0. If so, it means that the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point are correct, and the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point are recorded as the target initial value of the permittivity ε0 t0 and the target conductivity initial value σ t0 , without participating in the iteration, the iterative equation can be constrained by changing the unknown quantity to the known quantity K. The constraint term K is a vector, which is composed of the target initial value of the permittivity ε0 and the initial value of the conductivity σ0. The specific composition method is calculated by the equation deformation decomposition principle of the partial differential equation and the matrix equation, which will not be repeated here. If not, it means that the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point are unreasonable and need to participate in the iteration. The initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point are recorded as the iterative initial value of the permittivity ε0 d0 and the initial value of conductivity σ d0 , is substituted into the iterative equation as the quantity to be solved. The initial value of the iterative permittivity ε involved in the iteration is d0 and the initial value of conductivity σ d0 The scale is n x *n y *n z , n x 、ny 、n z Represents the number of pixels in the x, y, and z directions, so the number of pixels in the entire space is n x *n y *n z .
[0151] In the above step S34, the number of target initial conductivity values is the same as the number of target initial permittivity values, so the constant p can also be the number of target initial conductivity values. For example, 1000 reasonable target initial permittivity values and target initial conductivity values are screened out in steps S31 and S32, that is, the number of target initial permittivity values and target initial conductivity values is both 1000, so the constant p is 1000.
[0152] In the above step S4, based on the launch site and the point data information of the launch site, the high-order nonlinear nonhomogeneous partial differential equation derived from the original Maxwell equations and the Helmholtz equation is instantaneously linearized by the improved Newton iteration method to obtain the iterative formula, and the constraints K, R1 and R2 are added according to the physical model and the initial value of the permittivity ε0 and the initial value of the conductivity σ0 in reality to obtain the iterative equation
[0153] Furthermore, in the iterative equation, f(u) is as shown in formula (2):
[0154] f(u)=Au+Cu 2 (2)
[0155] Where u 2 is a column vector, u 2 The general formula for each term is A is a matrix determined by the launch site and the point data information of the launch site, as shown in formula (3):
[0156]
[0157] In the formula, the variable As shown in formula (4), the variable γ is shown in formula (5):
[0158]
[0159]
[0160] C is a matrix determined by the launch site and the point data information of the launch site, as shown in formula (6):
[0161]
[0162] Furthermore, in the iterative equation, the damping factor vector α is calculated as follows:
[0163] Build an n x * y * z The matrix Q of , the value of each point of Q is shown in formula (7):
[0164]
[0165] Where σ m and ε m are all parameters, (x i ,y i ,z i ) corresponds one-to-one with the subscript i, where i represents the i-th point;
[0166] Suppose there is a matrix V. For each point v in V, there is a corresponding m*m matrix T in Q with the coordinate of v as the center. The value of v is shown in formula (23):
[0167] v=min(T) (8)
[0168] The damping factor vector α is obtained by mapping the matrix V one by one to the damping factor vector α.
[0169] The element v of the matrix V corresponds to each point in space. The matrix T is equivalent to a scanning window. Its number of rows and columns m can be determined according to actual conditions. Preferably, m is 5. For each element v in the 5*5 matrix T, a 5*5 matrix T is constructed with the corresponding v coordinate point in the matrix Q as the center point, the center point and two points in each direction of the center point. The minimum value in the constructed 5*5 matrix T is found as the value of the element v. In this way, the value of each element v is continuously obtained, and then the matrix V is obtained. And from the above process, it can be seen that the matrix V will change at each iteration. Therefore, the matrix V needs to be recalculated after each iteration to recalculate the damping factor vector α.
[0170] In step S5 above, a reasonable initial value is used in the iterative method to accelerate convergence. Therefore, by determining the initial value of the iterative equation based on the initial values of the permittivity and conductivity at each location in space, the convergence rate can be improved. However, the magnetic resonance device used in the present invention is an orthogonal coil, so the component of the transmission field B in the z direction can be ignored. Therefore, the original equation involved is shown in formula (23):
[0171]
[0172] Where i represents an imaginary number, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, μ0 is the magnetic permeability, and B′ x is the partial derivative of the emission field B in the yz plane along the x-axis, B′ y is the partial derivative of the emission field B in the xz plane along the y axis, B′ zis the partial derivative of the emission field B in the xy plane along the z axis;
[0173] According to the Helmholtz EPT method, it is assumed that the electrical characteristic parameters are uniformly distributed in space, that is, the derivatives of the electrical characteristic parameters at each point in space are all 0 along all directions, so that formula (23) can be simplified as shown in formula (24):
[0174]
[0175] Formula (25) can be derived from formula (24):
[0176]
[0177] Then, the initial value u0=σ0+iωε0 required for the iterative equation is obtained. When constructing the initial value u0, it is transformed from a matrix composed of the initial values of the permittivity and the initial values of the conductivity at each point in the space into a vector u0, which has a one-to-one correspondence. The initial value u0 is substituted into the iterative equation for iterative calculation. After several iterations, the permittivity and conductivity at each point in the space are obtained.
[0178] Furthermore, based on the initial values of the permittivity and the conductivity at each point in the space, the initial values of the iterative equation are determined, the initial values are substituted into the iterative equation for iterative calculation, and the permittivity and conductivity at each point in the space are obtained after several iterations. The steps include:
[0179] S51. Preprocess the initial value of the iterative permittivity and the initial value of the iterative conductivity. Based on the preprocessed initial value of the iterative permittivity and the initial value of the iterative conductivity, calculate the initial value u0 of the iterative equation according to formula (15):
[0180] u0=σ0+iωε0 (15)
[0181] Where, σ0 is the initial value of iterative conductivity, ε0 is the initial value of iterative permittivity;
[0182] S52, substituting the initial value u0 into the iterative equation to obtain a new initial value, and determining whether the average value of each item in the new initial value is less than a preset iteration threshold;
[0183] S53, if not, repeat step S52;
[0184] S54. If yes, the permittivity and conductivity at each point in the space are obtained according to the new initial value after iteration, the target permittivity initial value, and the target conductivity initial value.
[0185] In the above steps S51 to S54, the initial value of the iterative permittivity ε that needs to be iterated is screened out by steps S31 to S33. d0 and the initial value of conductivity σ d0, according to formula (15), the initial value u0 is formed, and the initial value u0 is substituted into the iterative equation to obtain a new initial value u1. Whether the average value of each item in the new initial value u1 is less than the preset iteration threshold, it is judged whether the mean value of f(u1) is less than the iteration threshold. If not, the new initial value u1 is substituted into the iterative equation to obtain a new initial value u2. This is repeated until the initial value u of a certain generation is i If the average value of each item in is less than the preset iteration threshold, the iteration ends and u is obtained. i Since the initial value u0 of the component is transformed from a matrix to a vector, there is a one-to-one correspondence, so the initial value of each iterative permittivity ε can be recorded in advance d0 and the initial value of conductivity σ d0 The mapping information and the corresponding coordinates in space are also recorded in advance, and the target initial value of the permittivity ε that is not required to participate in the iteration is screened out from step S31 to step S33. t0 and the target conductivity initial value σ t0 The mapping information of u is calculated. i After that, based on u i , initial value of target permittivity ε t0 and the target conductivity initial value σ t0 , according to the coordinate information recorded in advance, it is restored into a matrix form, and the permittivity value ε and conductivity value σ at each point in the space can be obtained.
[0186] In the above step S6, the permittivity value ε and the conductivity value σ at each point in the space are the information of each point on the image, and an electrical property tomographic image can be output.
[0187] The following is a comparison of the results of the instantaneous linear electrical property tomography method of the present invention (hereinafter referred to as the present invention) with the Helmholtz EPT and the improved crEPT method:
[0188] The reference image is a Duke human brain tumor simulation model. Peak signal-to-noise ratio (PSNR) and structural similarity were used as evaluation metrics to compare the results of the present invention with those of the Helmholtz EPT and modified crEPT methods. Higher PSNR and structural similarity values indicate closer alignment of the reconstructed electrical tomography images with the reference image. Table 1 shows the quantitative metrics for the results obtained by the present invention, Helmholtz EPT, and modified crEPT methods.
[0189] Table 1
[0190] method name Peak signal-to-noise ratio Structural similarity Helmholtz EPT u 8.6151 0.6825 Modified crEPT u 15.5295 0.8178 The present invention u 18.2423 0.8684
[0191] As can be seen from Table 1, in this example, the present invention is superior to the Helmholtz EPT and improved crEPT methods in terms of both peak signal-to-noise ratio and structural similarity.
[0192] In addition, please refer to Figure 2 and Figure 3 , Figure 2 Comparison of the reconstructed conductivity (top) and permittivity (bottom) images of the Duke human brain tumor simulation model using the instantaneous linearization-based electrical property tomography method of the present invention, the Helmholtz EPT method, and the modified crEPT method. Figure 3 The results of reconstructing conductivity (top) and permittivity (bottom) images of another Duke human brain tumor simulation model using the improved crEPT method based on the instantaneous linearized electrical property tomography method and the Helmholtz EPT method of the present invention. Figure 2 and Figure 3 It can be seen that the image of the present invention is visually closer to the reference image. Therefore, the effect of the electrical property tomography image obtained by the present invention is significantly better than that of the Helmholtz EPT and the improved crEPT method.
[0193] Please refer to Figure 4 , Figure 4 The present invention also provides a module schematic diagram of an electrical property tomography system based on instantaneous linearization. The present invention also provides an electrical property tomography system based on instantaneous linearization, comprising:
[0194] The measurement module 1 is used to measure the magnetic resonance radio frequency field and calculate the transmission field based on the radio frequency field;
[0195] A first calculation module 2 is configured to calculate the point data information of the launch field using a numerical method, where the point data information includes the first-order derivative and the second-order derivative of the launch field along the x-axis, y-axis, and z-axis at each point in space;
[0196] The second calculation module 3 is used to calculate the initial value of the permittivity and the initial value of the conductivity at each point in space based on the launch site and the point data information of the launch site according to Maxwell's equations, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space;
[0197] Module 4 is obtained, which is used to obtain an iterative equation based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1):
[0198]
[0199] Where u n represents u of the nth generation, i.e., the iterative value obtained after the nth iteration, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector;
[0200] Iterative module 5 is used to determine the initial values of the iterative equation based on the initial values of the permittivity and the conductivity at each point in the space, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations;
[0201] The output module 6 is used to output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
[0202] In the measurement module 1, a nuclear magnetic resonance device, an orthogonal coil, and a specific B1 scanning sequence are used to obtain the radio frequency field B in the orthogonal coil measurement space, that is, the magnetic resonance radio frequency field B. The transmit field B is then calculated based on the obtained radio frequency field B.
[0203] Furthermore, the measurement module 1 includes:
[0204] The measurement submodule measures the components B of the magnetic resonance radio frequency field along the x-axis, y-axis and z-axis directions of the spatial coordinate system. x 、B y 、B z ;
[0205] The first calculation submodule is used to calculate the transmission field B according to formula (22):
[0206]
[0207] Where i is a complex unit.
[0208] In the above measurement submodule to the first calculation submodule, the unique vector of the radio frequency field B measured by the instrument is (B x ,B y ,B z ), where B x 、B y 、B z are the components of the radio frequency field B along the x-axis, y-axis and z-axis of the spatial coordinate system. The transmitting field B is given by the formula calculate.
[0209] The numerical method can adopt the finite difference method. Based on the launch field B, the first calculation module 2 calculates the first-order derivative and the second-order derivative of the launch field along the x-axis, y-axis and z-axis at each point in space through the finite difference method.
[0210] Furthermore, the first calculation module 2 includes:
[0211] The second calculation submodule is used to calculate the first-order derivative of the emission field along the x-axis at each point in space according to formula (16):
[0212]
[0213] The third calculation submodule is used to calculate the first-order derivative of the emission field along the y-axis at each point in space according to formula (17):
[0214]
[0215] The fourth calculation submodule is used to calculate the first-order derivative of the launch field along the z-axis at each point in space according to formula (18):
[0216]
[0217] The fifth calculation submodule is used to calculate the second-order derivative of the emission field along the x-axis at each point in space according to formula (19):
[0218]
[0219] The sixth calculation submodule is used to calculate the second-order derivative of the launch field along the y-axis at each point in space according to formula (20):
[0220]
[0221] The seventh calculation submodule is used to calculate the second-order derivative of the launch field along the z-axis at each point in space according to formula (21):
[0222]
[0223] Where dx, dy, and dz are the lengths of the unit grid in the x-axis, y-axis, and z-axis directions after spatial discretization, respectively, and i is a subscript indicating a specific point in space.
[0224] The second calculation module 3 calculates the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at each point in the space based on the launch field and the point data information of the launch field and according to the principle of MREPT (magnetic resonance electrical properties, magnetic resonance tomography of electrical properties of human tissue), and then processes the calculated initial value of the permittivity ε0 and the initial value of the conductivity σ0 at each point in the space according to the preset processing rules to determine the constraints K, R1 and R2.
[0225] Furthermore, the second calculation module 3 includes:
[0226] A judgment submodule is used to judge whether the initial values of the permittivity and the initial values of the conductivity at each point in the space are 0 along all directions;
[0227] a first marking submodule, configured to, if the judgment of the judgment submodule is yes, record the initial permittivity value and the initial conductivity value as the target permittivity value and the target conductivity value, exclude the target permittivity value and the target conductivity value from the iteration, and construct a constraint term K based on the target permittivity value and the target conductivity value;
[0228] The second marking submodule is used to record the initial value of the permittivity and the initial value of the conductivity as the iterative initial value of the permittivity and the iterative initial value of the conductivity if the judgment of the judgment submodule is no;
[0229] The first constraint module is used to constrain the item R1 to be of size 1*(n x *n y *n z ), the general formula of each item in the constraint term R1 is shown in formula (9):
[0230]
[0231] Where λ R1 is the coefficient of constraint R1, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, and n0 and (x0, y0, z0) satisfy the relationship shown in formula (10):
[0232] n0=n z *n y *x0+n z *y0+z0 (10)
[0233] f R1 As shown in formula (11):
[0234]
[0235] Where θ r It is a parameter, and its value is usually 0.2;
[0236] Constant term θ t0 As shown in formula (12):
[0237]
[0238] Where, the constant p is the number of target permittivity initial values;
[0239] The first constraint module is used to constrain the item R2 to be a size of 1*(n x *n y *n z ), the constraint term R2 is shown in formula (13):
[0240]
[0241] Where λ R2 is the coefficient of constraint R2, is the Laplace operator, u is a scale of 1*(n x *n y *n z ), each term of u is expressed as follows:
[0242]
[0243] Referring to the Helmholtz EPT method, it is assumed that the electrical characteristic parameters are uniformly distributed in space, that is, the derivatives of the electrical characteristic parameters along all directions are 0 at every point in space. This allows the partial differential equation to be converted into a linear equation, facilitating the acquisition of initial values for iteration, where the electrical characteristic parameters include permittivity and conductivity. However, in the real spatial distribution, some initial permittivity values ε0 and conductivity values σ0 do satisfy the assumption that the derivatives along all directions are equal to 0. These initial permittivity values ε0 and conductivity values σ0 are already very accurate and do not need to be solved again. Instead, they can be converted from unsolved quantities to known quantities to constrain the iterative equation itself, thus reducing the computational complexity and, therefore, the amount of computation required for subsequent iterative calculations. Therefore, based on the assumptions of the Helmholtz EPT method, the judgment submodule analyzes whether the initial values of the permittivity ε0 and the initial values of the conductivity σ0 at each point in space are 0 along the diagonals of all directions. If so, it means that the initial values of the permittivity ε0 and the initial values of the conductivity σ0 at the point are correct. The first marking module records the initial values of the permittivity ε0 and the initial values of the conductivity σ0 at the point as the target initial values of the permittivity ε0 t0 and the target conductivity initial value σ t0 , without participating in the iteration, the iterative equation can be constrained by changing the quantity to be determined to the known quantity K. The constraint term K is a vector, which is composed of the target initial value of the permittivity ε0 and the initial value of the conductivity σ0. The specific composition method is calculated by the equation deformation decomposition principle of the partial differential equation and the matrix equation, which will not be repeated here. If not, it means that the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point are unreasonable and need to participate in the iteration. The second marking module records the initial value of the permittivity ε0 and the initial value of the conductivity σ0 at this point as the iterative initial value of the permittivity ε0 d0 and the initial value of conductivity σ d0 , is substituted into the iterative equation as the quantity to be solved. The initial value of the iterative permittivity ε involved in the iteration is d0 and the initial value of conductivity σ d0 The scale is n x *n y *n z , n x 、n y 、n zRepresents the number of pixels in the x, y, and z directions, so the number of pixels in the entire space is n x *n y *n z .
[0244] In the first constraint module, the number of target initial conductivity values is the same as the number of target initial permittivity values. Therefore, the constant p can also be the number of target initial conductivity values. For example, the judgment submodule and the first marking submodule screen out 1000 reasonable target initial permittivity values and target initial conductivity values. That is, the number of target initial permittivity values and target initial conductivity values is 1000, so the constant p is 1000.
[0245] Module 4 obtains the iterative formula by instantaneously linearizing the high-order nonlinear nonhomogeneous partial differential equation derived from the original Maxwell equations and the Helmholtz equation based on the launch site and the point data information of the launch site, and adds the constraints K, R1 and R2 according to the physical model and the initial value of the permittivity ε0 and the initial value of the conductivity σ0 in reality, and obtains the iterative equation
[0246] Furthermore, in the iterative equation, f(u) is as shown in formula (2):
[0247] f(u)=Au+Cu 2 (2)
[0248] Where u 2 is a column vector, u 2 The general formula for each term is A is a matrix determined by the launch site and the point data information of the launch site, as shown in formula (3):
[0249]
[0250] In the formula, the variable As shown in formula (4), the variable γ is shown in formula (5):
[0251]
[0252]
[0253] C is a matrix determined by the launch site and the point data information of the launch site, as shown in formula (6):
[0254]
[0255] Furthermore, in the iterative equation, the damping factor vector α is calculated as follows:
[0256] Build an n x *y * z The matrix Q of , the value of each point of Q is shown in formula (7):
[0257]
[0258] Where σ m and ε m are all parameters, (x i ,y i ,z i ) corresponds one-to-one with the subscript i, where i represents the i-th point;
[0259] Suppose there is a matrix V. For each point v in V, there is a corresponding m*m matrix T in Q with the coordinate of v as the center. The value of v is shown in formula (23):
[0260] v=min(T) (8)
[0261] The damping factor vector α is obtained by mapping the matrix V one by one to the damping factor vector α.
[0262] The element v of the matrix V corresponds to each point in space. The matrix T is equivalent to a scanning window. Its number of rows and columns m can be determined according to actual conditions. Preferably, m is 5. For each element v in the 5*5 matrix T, a 5*5 matrix T is constructed with the corresponding v coordinate point in the matrix Q as the center point, the center point and two points in each direction of the center point. The minimum value in the constructed 5*5 matrix T is found as the value of the element v. In this way, the value of each element v is continuously obtained, and then the matrix V is obtained. And from the above process, it can be seen that the matrix V will change at each iteration. Therefore, the matrix V needs to be recalculated after each iteration to recalculate the damping factor vector α.
[0263] In the above-mentioned iterative module 5, a reasonable initial value is used in the iterative method to accelerate convergence. Therefore, by determining the initial value of the iterative equation based on the initial values of the permittivity and conductivity at each location in space, the convergence speed can be improved. However, the magnetic resonance device used in the present invention is an orthogonal coil, so the component of the transmission field B in the z direction can be ignored. Therefore, the original equation involved is shown in formula (23):
[0264]
[0265] Where i represents an imaginary number, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, μ0 is the magnetic permeability, and B′ x is the partial derivative of the emission field B in the yz plane along the x-axis, B′ y is the partial derivative of the emission field B in the xz plane along the y axis, B′ z is the partial derivative of the emission field B in the xy plane along the z axis;
[0266] According to the Helmholtz EPT method, it is assumed that the electrical characteristic parameters are uniformly distributed in space, that is, the derivatives of the electrical characteristic parameters at each point in space are all 0 along all directions, so that formula (23) can be simplified as shown in formula (24):
[0267]
[0268] Formula (25) can be derived from formula (24):
[0269]
[0270] Then, the initial value u0=σ0+iωε0 required for the iterative equation is obtained. When constructing the initial value u0, it is transformed from a matrix composed of the initial values of the permittivity and the initial values of the conductivity at each point in the space into a vector u0, which has a one-to-one correspondence. The initial value u0 is substituted into the iterative equation for iterative calculation. After several iterations, the permittivity and conductivity at each point in the space are obtained.
[0271] Furthermore, the iteration module 5 includes:
[0272] The processing submodule is used to preprocess the initial value of the iterative permittivity and the initial value of the iterative conductivity, and calculate the initial value u0 of the iterative equation according to formula (15) based on the preprocessed initial value of the iterative permittivity and the initial value of the iterative conductivity:
[0273] u0=σ0+iωε0 (15)
[0274] Where, σ0 is the initial value of iterative conductivity, ε0 is the initial value of iterative permittivity;
[0275] The substitution module is used to substitute the initial value u0 into the iterative equation to obtain a new initial value and determine whether the average value of each item in the new initial value is less than a preset iteration threshold;
[0276] Repeat submodule, for repeating step S52 if no;
[0277] A submodule is obtained for obtaining the permittivity and conductivity at each point in the space according to the new initial value after iteration, the target permittivity initial value, and the target conductivity initial value.
[0278] In the above processing submodule to the obtaining submodule, the initial value ε of the iterative capacitance ratio that needs to be iterated is screened out by the judging submodule, the first marking submodule and the second marking submodule. d0 and the initial value of conductivity σ d0, according to formula (15), the initial value u0 is formed, and the initial value u0 is substituted into the iterative equation to obtain a new initial value u1. Whether the average value of each item in the new initial value u1 is less than the preset iteration threshold, it is judged whether the mean value of f(u1) is less than the iteration threshold. If not, the new initial value u1 is substituted into the iterative equation to obtain a new initial value u2. This is repeated until the initial value u of a certain generation is i If the average value of each item in is less than the preset iteration threshold, the iteration ends and u is obtained. i Since the initial value u0 of the component is transformed from a matrix to a vector, there is a one-to-one correspondence, so the initial value of each iterative permittivity ε can be recorded in advance d0 and the initial value of conductivity σ d0 The mapping information and the corresponding coordinates in space are also recorded in advance, and the target initial value of the permittivity ε that is not required to participate in the iteration is screened out from step S31 to step S33. t0 and the target conductivity initial value σ t0 The mapping information of u is calculated. i After that, based on u i , initial value of target permittivity ε t0 and the target conductivity initial value σ t0 , according to the coordinate information recorded in advance, it is restored into a matrix form, and the permittivity value ε and conductivity value σ at each point in the space can be obtained.
[0279] In the output module 6 , the permittivity value ε and the conductivity value σ at each point in the space are the information of each point on the image, and an electrical property tomographic image can be output.
[0280] Please refer to Figure 5 , Figure 5 This is a schematic structural block diagram of an electronic device according to an embodiment of the present invention. An embodiment of the present invention further proposes an electronic device 1001, comprising a memory 1003 and a processor 1002, wherein the memory 1003 stores a computer program 1004, and when the processor 1002 executes the computer program 1004, the steps of any of the above-mentioned electrical property tomography methods based on instantaneous linearization are implemented, including: S1, measuring the magnetic resonance radio frequency field, and calculating the transmission field based on the radio frequency field; S2, using a numerical method to calculate the point data information of the transmission field, the point data information including the first-order derivative and the second-order derivative of the transmission field along the x-axis, y-axis and z-axis directions at each point in space; S3, based on the transmission field and the point data information of the transmission field, according to Maxwell's equations, calculating the initial value of the permittivity and the initial value of the conductivity at each point in space, and determining the constraints K, R1 and R2 according to the initial value of the permittivity and the initial value of the conductivity at each point in space; S4, based on the transmission field, the point data information of the transmission field and the constraints K, R1 and R2, according to Maxwell's equations combined with Newton's iteration method, obtaining the iterative equation S5. Based on the initial values of the permittivity and conductivity at each point in the space, determine the initial values of the iterative equation, substitute the initial values into the iterative equation for iterative calculation, and after several iterations, obtain the permittivity and conductivity at each point in the space; S6. Output the electrical characteristic tomographic image based on the permittivity and conductivity at each point in the space.
[0281] Please refer to Figure 6 , Figure 6 The present invention is a block diagram of a computer-readable storage medium according to an embodiment of the present invention. An embodiment of the present invention further provides a computer-readable storage medium 2001, on which a computer program 1004 is stored. When the computer program 1004 is executed by the processor 1002, the computer program 1004 implements the steps of any of the above-mentioned electrical characteristic tomography methods based on instantaneous linearization, including: S1, measuring the magnetic resonance radio frequency field, and calculating the transmission field based on the radio frequency field; S2, using a numerical method to calculate the point data information of the transmission field, the point data information including the first-order derivative and the second-order derivative of the transmission field along the x-axis, y-axis and z-axis at each point in space; S3, based on the transmission field and the point data information of the transmission field, according to Maxwell's equations, calculate the initial value of the permittivity and the initial value of the conductivity at each point in space, and determine the constraints K, R1 and R2 according to the initial value of the permittivity and the initial value of the conductivity at each point in space; S4, based on the transmission field, the point data information of the transmission field and the constraints K, R1 and R2, according to Maxwell's equations combined with Newton's iteration method, obtain the iterative equation S5. Based on the initial values of the permittivity and conductivity at each point in the space, determine the initial values of the iterative equation, substitute the initial values into the iterative equation for iterative calculation, and after several iterations, obtain the permittivity and conductivity at each point in the space; S6. Output the electrical characteristic tomographic image based on the permittivity and conductivity at each point in the space.
[0282] Compared with the existing technology, the beneficial effects of the present invention are: the high-order nonlinear nonhomogeneous partial differential equations derived from the original Maxwell equations and the Helmholtz equations are instantaneously linearized through the improved Newton iteration method to obtain iterative equations, and constraint terms are added according to the physical model and actual conditions to form a complete iterative equation with initial values, thereby obtaining an electrical characteristic distribution solution with higher accuracy and more in line with the actual situation, thereby being able to obtain electrical characteristic tomography more accurately while minimizing the simplification or change of the problem model.
[0283] Those skilled in the art will understand that all or part of the processes in the above-mentioned embodiment methods can be implemented by instructing the relevant hardware through a computer program, and the computer program can be stored in a non-volatile computer-readable storage medium. When the computer program is executed, it can include the processes of the embodiments of the above-mentioned methods. Among them, any reference to memory, storage, database or other media provided in this application and used in the embodiments may include non-volatile and / or volatile memory. Non-volatile memory may include read-only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM) or flash memory. Volatile memory may include random access memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in many forms, such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), double data rate SDRAM (SSRSDRAM), enhanced SDRAM (ESDRAM), Synchronous Link DRAM (SLDRAM), Rambus direct RAM (RDRAM), direct RAM bus dynamic RAM (DRDRAM), and RAM bus dynamic RAM (RDRAM), etc.
[0284] It should be noted that, in this document, the terms "comprises," "includes," or any other variations thereof are intended to encompass non-exclusive inclusion, such that a process, apparatus, article, or method comprising a series of elements includes not only those elements but also other elements not explicitly listed, or elements inherent to such process, apparatus, article, or method. In the absence of further limitations, an element defined by the phrase "comprising a ..." does not exclude the presence of other identical elements in the process, apparatus, article, or method comprising the element.
[0285] The present invention is not limited to the above-mentioned embodiments. If various changes or modifications of the present invention do not depart from the spirit and scope of the present invention, and if these changes and modifications fall within the scope of the claims of the present invention and equivalent technologies, the present invention is also intended to include these changes and modifications.
Claims
1. An electrical property tomography method based on instantaneous linearization, characterized in that: The following steps are involved: S1. Measure the magnetic resonance radio frequency field and calculate the transmit field based on the radio frequency field; S2. Calculate point data information of the launch field using a numerical method, where the point data information includes first-order derivatives and second-order derivatives of the launch field along the x-axis, y-axis, and z-axis at each point in space; S3. Based on the launch site and the point data information of the launch site, according to Maxwell's equations, calculate the initial value of the permittivity and the initial value of the conductivity at each point in space, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space; S4. Based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, an iterative equation is obtained according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1): Where u n represents u of the nth generation, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector; S5. Based on the initial values of the permittivity and conductivity at each point in the space, determine the initial values of the iterative equation, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations; S6. Output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
2. The electrical property tomography method based on instantaneous linearization according to claim 1, characterized in that: In the iterative equation, f(u) is as shown in formula (2): f(u)=Au+Cu 2 (2) Where u 2 is a column vector, u 2 The general formula for each term is A is a matrix, as shown in formula (3): In the formula, the variable As shown in formula (4), the variable γ is shown in formula (5): C is a matrix, as shown in formula (6):
3. The electrical property tomography method based on instantaneous linearization according to claim 1, characterized in that: In the iterative equation, the damping factor vector α is calculated as follows: Build an n x *n y *n z The matrix Q of , the value of each point of Q is shown in formula (7): Where σ m and ε m All are parameters; Suppose there is a matrix V. For each point v in V, there is a corresponding m*m matrix T in Q with the coordinate of v as the center. The value of v is shown in formula (23): v=min(T) (8) The damping factor vector α is obtained by mapping the matrix V one by one to the damping factor vector α.
4. The electrical property tomography method based on instantaneous linearization according to claim 1, characterized in that: The step of determining the constraints K, R1 and R2 based on the initial values of the permittivity and conductivity at each point in the space includes: For the initial value of the permittivity and the initial value of the conductivity at each point in the space, determine whether the initial value of the permittivity and the initial value of the conductivity along all directions are 0; If so, the initial values of permittivity and conductivity are recorded as the target initial values of permittivity and conductivity, and the target initial values of permittivity and conductivity do not participate in the iteration. The constraint term K is constructed based on the target initial values of permittivity and conductivity. If not, the initial values of permittivity and conductivity are recorded as the initial values of iterative permittivity and iterative conductivity; The constraint R1 is of size 1*(n x *n y *n z ), the general formula of each item in the constraint term R1 is shown in formula (9): Where λ R1 is the coefficient of constraint R1, ω is the Larmor precession frequency corresponding to the nuclear magnetic resonance device, and n0 and (x0, y0, z0) satisfy the relationship shown in formula (10): n0=n z *n y *x0+n z *y0+z0 (10) f R1 As shown in formula (11): Where θ r is a parameter; Constant term θ t0 As shown in formula (12): Where, the constant p is the number of target permittivity initial values; The constraint R2 is a scale of 1*(n x *n y *n z ), the constraint term R2 is shown in formula (13): Where λ R2 is the coefficient of constraint R2, is the Laplace operator, u is a scale of 1*(n x *n y *n z ), each term of u is expressed as follows:
5. The electrical property tomography method based on instantaneous linearization according to claim 4, characterized in that: The steps of determining initial values of an iterative equation based on initial values of the permittivity and conductivity at each point in the space, substituting the initial values into the iterative equation for iterative calculation, and obtaining the permittivity and conductivity at each point in the space after several iterations include: S51. Preprocess the initial value of the iterative permittivity and the initial value of the iterative conductivity. Based on the preprocessed initial value of the iterative permittivity and the initial value of the iterative conductivity, calculate the initial value u0 of the iterative equation according to formula (15): u0=σ0+iωε0 (15) Where, σ0 is the initial value of iterative conductivity, ε0 is the initial value of iterative permittivity; S52, substituting the initial value u0 into the iterative equation to obtain a new initial value, and determining whether the average value of each item in the new initial value is less than a preset iteration threshold; S53, if not, repeat step S52; S54. If yes, the permittivity and conductivity at each point in the space are obtained according to the new initial value after iteration, the target permittivity initial value, and the target conductivity initial value.
6. The electrical property tomography method based on instantaneous linearization according to claim 1, characterized in that: The step of calculating the point data information of the launch site using a numerical method includes: According to formula (16), the first-order derivative of the launch field along the x-axis at each point in space is calculated as follows: According to formula (17), the first-order derivative of the launch field along the y-axis at each point in space is calculated as follows: According to formula (18), the first-order derivative of the launch field along the z-axis at each point in space is calculated as follows: According to formula (19), the second-order derivative of the launch field along the x-axis at each point in space is calculated as According to formula (20), the second-order derivative of the launch field along the y-axis at each point in space is calculated as According to formula (21), the second-order derivative of the launch field along the z-axis at each point in space is calculated as Among them, dx, dy, and dz are the lengths of the unit grid in the x-axis, y-axis, and z-axis directions after spatial discretization, respectively.
7. The electrical property tomography method based on instantaneous linearization according to claim 1, characterized in that: The step of measuring the magnetic resonance radio frequency field and calculating the transmission field based on the radio frequency field includes: Measure the components B of the magnetic resonance radio frequency field along the x-axis, y-axis, and z-axis of the spatial coordinate system x 、B y 、B z ; Calculate the launch field B according to formula (22): Where i is a complex unit.
8. An electrical property tomography system based on instantaneous linearization, characterized in that: include: A measurement module is used to measure the magnetic resonance radio frequency field and calculate the transmission field based on the radio frequency field; A first calculation module is configured to calculate point data information of the launch field using a numerical method, wherein the point data information includes first-order derivatives and second-order derivatives of the launch field along the x-axis, y-axis, and z-axis at each point in space; The second calculation module is used to calculate the initial value of the permittivity and the initial value of the conductivity at each point in space based on the launch site and the point data information of the launch site according to Maxwell's equations, and determine the constraints K, R1 and R2 based on the initial value of the permittivity and the initial value of the conductivity at each point in space; The module is used to obtain an iterative equation based on the launch site, the point data information of the launch site, and the constraints K, R1, and R2, according to the Maxwell equations combined with the Newton iteration method. The iterative equation is shown in formula (1): Where u n represents u of the nth generation, f(u) is a function determined by the point data information of the launch site B, and α is the damping factor vector; An iterative module is used to determine the initial values of the iterative equation based on the initial values of the permittivity and conductivity at each point in the space, substitute the initial values into the iterative equation for iterative calculation, and obtain the permittivity and conductivity at each point in the space after several iterations; The output module is used to output an electrical property tomographic image based on the permittivity and conductivity at each point in space.
9. An electronic device comprising a memory and a processor, wherein the memory stores a computer program, wherein: When the processor executes the computer program, the steps of the method according to any one of claims 1 to 7 are implemented.
10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the steps of the method according to any one of claims 1 to 7 are implemented.
Citation Information
Patent Citations
System and method for data reconstruction in soft-field tomography
CN103099616A
System and method for determining fluid and tissue volume estimates using electrical characteristic tomography
CN114270397A