CN105224764B - Bone modeling and simulation method - Google Patents
Bone modeling and simulation method Download PDFInfo
- Publication number
- CN105224764B CN105224764B CN201510689648.9A CN201510689648A CN105224764B CN 105224764 B CN105224764 B CN 105224764B CN 201510689648 A CN201510689648 A CN 201510689648A CN 105224764 B CN105224764 B CN 105224764B
- Authority
- CN
- China
- Prior art keywords
- point
- bone
- spline
- value
- parameter value
- 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.)
- Expired - Fee Related
Links
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Image Generation (AREA)
Abstract
The present invention provides a kind of bone modeling and simulation methods, comprising: step 1, carries out CT scan to bone to obtain image data, obtains point cloud data；Step 2, basal plane projection is carried out according to the profile of point cloud data, obtains the parameter value and granule surface contral point of point cloud data；Step 3, body interpolation is carried out using granule surface contral point as known conditions, obtains body parameterized model using discrete Coons body interpolation algorithm；Step 4, multiple sampled points are set, and body parameterized model is marked off into multiple nodes, the gray value of granule surface contral point and interior control point is acquired by principle of least square method and quadratic form optimization；Step 5, the skeleton model of bone is obtained according to the parameter value of sampled point and gray value；Step 6, it carries out grid subdivision and obtains unit grid, the element stiffness matrix of unit grid is attached in global stiffness matrix, the value of unknown quantity is then obtained, finally treated that result is shown to the value of unknown quantity.
Description
Technical field
The invention belongs to skeleton models to construct field, and in particular to one kind can construct more accurate bone modeling and imitate
True method.
Background technique
Emulation is using the essential process that occurs in model reproduction real system, with the development of science and technology, emulation technology
Using more and more extensive, the systems such as electrical, mechanical, chemical industry, waterpower, heating power are applied to, also include society, economic, biology doctor
The systems such as, management.In biomedical aspect, emulation has great importance on the medical history of the mankind.There is biological doctor
After learning emulation, doctor can reduce in this way the pain of patient, and cost by virtual sick body come lift technique
Low, repeatable, visual result.Human body collsion damage occurs most frequently.So having by the research to bone collsion damage
Very important medical value.Currently, carrying out l-G simulation test using by establishing skeleton model.
Generally, the model towards material reverse building product, reverse engineering includes that geometry reverse, material reverse and technique are anti-
It asks.For the material reverse of homogenous product, due to material property isotropism, it is easy to several main material properties
Parameter expresses the material property of entire product, and reverse focuses on the test of material property parameter, and bone belong to it is heterogeneous
Product, material property anisotropy, and the material property for being located at different spatial point may also be different, the material of bone
Performance and geometric position etc. have closely related.Traditional material detection means such as hardness measurement, metal lographic examination, chemical analysis, light
The methods of spectrum analysis is mainly for homogeneous material, and therefore, detection means was once becoming the difficult point of bone material reverse.In addition, material
The problem of expecting reverse is heterogeneous material design and representation, existing bone material expression are not too ideal.
Summary of the invention
The present invention is to carry out in order to solve the above problems, and it is an object of the present invention to provide a kind of can construct more accurate bone
Bone modeling and simulation method.
The present invention provides a kind of bone modeling and simulation methods, for carrying out modeling and simulation to human skeleton,
It is characterized in that, comprising the following steps: step 1, CT scan is carried out to obtain image data, then, to image data to bone
Point cloud data is obtained after carrying out image procossing；Step 2, according to the profile of point cloud data carry out basal plane projection, and with B-spline phase
It, then, will be in point cloud data to obtain the parameter value of point cloud data to that should be parameterized to the B-spline that point cloud data carries out
Scattered points carry out base curved surface projection and obtain the parameter value of the scattered points, then, scattered points are carried out B using principle of least square method
Spline surface is fitted and obtains granule surface contral point；Step 3, body interpolation is carried out using granule surface contral point as known conditions, utilized
Discrete Coons body interpolation algorithm obtains body parameterized model, and the granule surface contral point in body parameterized model corresponds to the surface of bone
Information, the interior control point in body parameterized model correspond to the material information of bone；Step 4, multiple sampled points are set, and by body
Parameterized model marks off multiple nodes, and then, searching and the nearest node of sampled point, the parameter value of the node is as sampled point
Initial parameter value, then near initial parameter value AMSAA discrete model region, it is nearest to find out distance sample in the zone
B-spline body node, using the B-spline body node as the parameter value of sampled point；Also, according to the parameter value of sampled point and pass through
Principle of least square method and quadratic form optimization acquire the gray value of granule surface contral point and interior control point；Step 5, according to sampled point
Parameter value and gray value obtain the skeleton model of bone；And step 6, grid subdivision is carried out according to B-spline and obtains unit grid,
Unknown quantity is expressed as form corresponding with the basic function of B-spline, then uses weak form and minimum potential energy principal will
The element stiffness matrix of unit grid is attached in global stiffness matrix, is then handled global stiffness matrix, is obtained unknown
The value of amount, finally to the value of unknown quantity, treated that result is shown.
It in bone modeling and simulation method provided by the invention, can also have the following features: wherein, step 1 packet
Containing the pretreatment for carrying out noise filtering, sharpening feature and edge enhancing to image data.
It in bone modeling and simulation method provided by the invention, can also have the following features: wherein, step 1 is also
It is set comprising the threshold values to image data, so that the contrast of the bone in image data increases, by image data
Interpolated value processing and define the range of threshold values, and principle is increased to being split in image data using region, thus
Obtain model sum of the grayscale values point cloud data.
In bone modeling and simulation method provided by the invention, it can also have the following features: wherein, in step 4,
Parameterizing body Model according to heterogeneous B-spline indicates are as follows:
In formula, Ni,p,Nj,q,Nk,rFor the B-spline basic function of 3D solid in three directions；H and λijkRespectively indicate section
The grayscale information of point and control point, and normalized has been carried out,
If sampled point is Ts={ xs,ys,zs,hs, i.e., each sampled point meets:
According to matrix properties, formula can also be write as following form:
The quantity of node is ten times of sampled point quantity.
In bone modeling and simulation method provided by the invention, it can also have the following features: wherein, in step 4,
The nearest node of distance sample is found using OBBs algorithm.
In bone modeling and simulation method provided by the invention, it can also have the following features: wherein, utilize minimum
Square law principle solves the optimization problem that following formula indicates:
N is set in formulaijks=NI, p(us)Nj,q(vs)Nk,r(ws), wherein Γ (λijk) it is to meet principle of least square method
λijk, it is rewritten into:
α=(i-1)+(j-1) × l+ (k-1) × l × m+1 is set, optimizing expression is obtained:
According to the property of B-spline basic function, obtain:
Also, last solution are as follows:
Due to λαConstraint more than or equal to 0 and less than or equal to 1, obtains control point gray value λijk。
In bone modeling and simulation method provided by the invention, it can also have the following features: wherein, in step 6,
Global stiffness matrix is handled, linear equation is solved:
Ku=F,
In formula, K is global stiffness matrix, and u is transposed matrix, and F is loading matrix, and wherein K is by assembling unit rigidity
Matrix KeIt obtains,
Element stiffness matrix is calculated with following formula:
In above formula:
Rs=Ni(ug)Nj(vg)Nk(wg),
In formula, S=i+j+k,
i∈[1,p+1],j∈[1,q+1],k∈[1,r+1]s∈[1,ne],ne=(p+1) (q+1) (r+1),
ug,vg,wgIt is calculated according to Gauss integration point selection rule common in finite element method,
Also, it sets
In formula, it is Poisson's ratio that E, which is Young's modulus, ν,.
The action and effect of invention
Related bone modeling and simulation method according to the present invention, by human skeleton carry out CT scan and into
Point cloud data is extracted after row pretreatment, then, basal plane projection is carried out according to the profile of point cloud data and obtains the ginseng of point cloud data
Numerical value, also, the scattered points in point cloud data are subjected to base curved surface projection and obtain the parameter value of scattered points, then, utilize minimum
Scattered points are carried out B-spline surface fitting and obtain granule surface contral point by square law principle；Utilize discrete Coons body interpolation algorithm
It obtains body parameterized model, finds the nearest node of sampled point, initial parameter value of the parameter value of the node as sampled point, so
Afterwards near initial parameter value AMSAA discrete model region, the B-spline body node nearest from sampled point is found out in the zone, by the B
Then parameter value of the batten body node as sampled point acquires granule surface contral by principle of least square method and quadratic form optimization
The gray value of point and interior control point, then obtains skeleton model, finally, carrying out grid subdivision according to B-spline obtains unit grid,
Unknown quantity is expressed as form corresponding with the basic function of B-spline, then using weak form and minimum potential energy principal
The element stiffness matrix of unit grid is attached in global stiffness matrix, then global stiffness matrix is handled, is obtained not
The value for the amount of knowing, finally to the value of unknown quantity, treated that result is shown, so, bone modeling and simulation method of the invention
More accurate skeleton model, and energy accurate expression real material distribution can be constructed, is more in line with reality to design
The bone product of situation.
Detailed description of the invention
Fig. 1 is the flow chart of bone modeling and simulation method in the embodiment of the present invention；
Fig. 2 is the volume mesh display figure of bone modeling and simulation method in the embodiment of the present invention；
Fig. 3 is the exterior view of the skeleton model of bone modeling and simulation method in the embodiment of the present invention；And
Fig. 4 is the strain analysis figure of bone modeling and simulation method in the embodiment of the present invention.
Specific embodiment
It is real below in order to be easy to understand the technical means, the creative features, the aims and the efficiencies achieved by the present invention
Example combination attached drawing is applied to be specifically addressed bone modeling and simulation method of the invention.
Fig. 1 is the flow chart of bone modeling and simulation method in the embodiment of the present invention.
As shown in Figure 1, in the present embodiment, bone modeling and simulation method passes through body parametric method implementation model geometry
Space construction realizes material space construction by the parametrization to gray value；By solving a quadratic form optimization problem, by material
Material space is embedded in geometric space, realizes the coupling of geometric space and material space, by waiting geometrical analysis to obtain the property of bone
Energy information, to be finally completed the building of body parametrization bone product model.
The present embodiment constructs the body parameterized model of bone using the thought of material reverse.Using CT scan data as several
The source of what and material information, constructs geometric space using the method for body parametrization, constructs material using the grayscale information of extraction
Space.Based on geometric space, material space is embedded into geometric space, realizes the coupling of geometric space and material space
Conjunction is emulated to obtain performance information finally by equal geometrical analysis, so that the moulding of heterogeneous product bone is realized, thus
Realize the moulding of true bone.
The data that the present embodiment is provided using CT scan bone are as the input data of true skeleton model, CT scan data
Geometry and material information can be exported simultaneously in the form of cloud.For geological information, body is rebuild using point cloud data and parameterizes mould
Type realizes universe body parametrization.Using only material information, it is able to achieve the 2.5 dimensional parameterizations expression of material space.Utilize material
Material space embed geometric space method, complete the material information body Parameter Expression of heterogeneous product, thus realize geometry and
The entity Parameter Expression of material.The accurate mould of skeleton that this method can not only provide parametrization, meeting actual distribution
Type, and model uses body Parameter Expression, provides effective representation aids for material reverse, and show and wait for subsequent body
Reliable model basis is established in geometrical analysis.Bone body parameterized model can not only accurately express the distribution of product real material,
Material redesign can be more easily carried out simultaneously, as changed the material distribution at control point, therefore this method can be designed and more be accorded with
Close the bone product of actual state.
Bone modeling and simulation method to human skeleton carry out modeling and simulation, bone modeling and simulation method it is specific
Steps are as follows:
Step S1 carries out CT scan to bone to obtain image data, then, mentions after pre-processing to image data
Get model sum of the grayscale values point cloud data；By then being pre-processed to image data to true bone progress CT scan,
In image pre-processing phase, needs to carry out picture noise filtering, sharpens the processing such as feature and edge enhancing, then, will own
The tomography picture that processing is completed is imported into sequence in Medical Image Processing software Mimics, passes through the threshold values to image data
It is set, so that bone portion and other contrast in tissue increase, convenient for accurately extracting skeleton model, for local threshold values
Very close to the method using manual identified.It is handled by interpolated value, defines the range of threshold values, principle is increased to not using region
It is split with region.After marking all surface profile points, processing to complete simultaneously, model gray value and point cloud number are extracted
According to subsequently into step S2.
Step S2, the surface fitting based on point cloud data, firstly, constructing a spline surface according to the profile of point cloud data
It is complete to put the projective parameter value on basal plane as the parameter value of data point by point cloud data to the basal plane projection as basal plane
It is parameterized at the B-spline of cloud；Then, construct B-spline surface parametric equation with the point cloud data that has parameterized, with etc. partial nodes
Basal plane is divided into different zones block, by the scattered points in point cloud data in region to base curved surface projection, calculates separately scattered points
At a distance from each equal partial nodes, the corresponding node of minimum range is the scattered points initial projections parameter, calculates the parameter of subpoint
It is worth the parameter value as scattered points.Then B-spline surface fitting, six clouds are carried out to scattered points using principle of least square method
The surface fitting of dough sheet obtains granule surface contral point after completing, and by increasing control points and number of nodes, realizes to surface fitting
The control of error, subsequently into step S3.
Step S3 is carried out body interpolation for granule surface contral point as known conditions, is obtained using discrete Coons body interpolation algorithm
Body parameterized model, the granule surface contral point of body parameterized model correspond to the surface information of bone, the interior control of body parameterized model
The corresponding bone of point and material information, subsequently into step S4.
Step S4, geometry and material Coupling method:
Homogeneous B-spline parametrization body Model is expressed as shown in formula (1)
In formula (1), Ni,p,Nj,q,Nk,rFor the B-spline basic function of 3D solid in three directions；H and λijkIt respectively indicates
Node T (u, v, w) and control point PijkThe grayscale information of (x, y, z), and carried out normalized.Each control point is come
It says, geometric coordinate (x, y, z) is known, gray value λijkIt is unknown and need solve, therefore solve λijkIt is geometry
Where the purpose and difficult point of material Coupling method, need according to the grayscale information of sampled point come the grayscale information at reverse control point.
If sampled point is Ts={ xs,ys,zs,hs, for each sampled point on slice, gray value hsIt is known
, i.e., each sampled point meets formula
According to matrix properties, formula can also be write as following form:
Corresponding parameter (the u of each sampled points,vs,ws) do not directly give, it needs to acquire by certain step.Cause
This, following steps can be divided by seeking step:
Step S4-1 seeks the parameter value (u of sampled points,vs,ws): multiple sampled points are set, and body is parameterized into mould
Type marks off multiple nodes, then, finds and the nearest node of sampled point, initial ginseng of the parameter value of the node as sampled point
Numerical value, then the further discrete local B-spline domain near initial parameter value, finds out the B nearest from sampled point in B-spline domain
Batten body node, using the B-spline body node as the parameter value of sampled point.
To guarantee to obtain preferable discrete parameter as a result, if the number of model sampled point is k, by whole individual parameter
Change model partition, which goes out 10*k node, to be advisable.For certain sampled point T on slices, find the corresponding parameter value of nearest node
(us,vs,ws) it is used as TsInitial parameter value in B-spline body.In the process, OBBs algorithm can be taken to accelerate operation speed
Degree, i.e., construct the bounding box of each infinitesimal ferritic first.If sampled point is fallen in some bounding box, the sampled point is being searched for
When nearest node, closest approach need to be only searched in the neighbouring bounding box of the bounding box.If sampled point fall in bounding box it
Outside, then the bounding box of search minimum distance is needed to search nearest point.Then in initial parameter value (us,vs,ws) in environs
Discrete local B-spline body domain again, is then searched for from this domain of individuals from sampled point T againsNearest B-spline body node.It will most
Whole parameter value (us,vs,ws) it is used as TsParameter value in B-spline body.Same method obtains being sliced upper each sampled point
Parameter value, subsequently into step S4-2.
Step S4-2 seeks control gray value λijk: due to meeting formula (3) at each sampled point, i.e. the equation left side can
It is directly obtained by sampling point information, equation the right (u because known tos,vs,ws) and the value of basic function can be calculated.If control point number is
M=l × m × n, sampled point number are N, and General N will be far longer than M.Using principle of least square method, which translates into one
The optimization problem that a optimization problem, i.e. solution following formula (4) indicate.
s.t.0≤λijk≤1 (4)
If Nijks=Ni,p(us)Nj,q(vs)Nk,r(ws), wherein Γ (λijk) and following Γ (λα) it is to meet least square method
The λ of principleijkAnd λαValue.
Formula (4) can be further rewritten as:
If remembering α=(i-1)+(j-1) × l+ (k-1) × l × m+1, and expansion (6), then further obtain shaped like (7) formula
Shown in optimizing expression:
The problem is converted into a quadratic form optimization problem, and due to the property of B-spline basic function, the solution of the system is centainly deposited
?.Formula (7) is further written as follow form:
Therefore, the last solution of formula (8) are as follows:
Due to λαThe presence of ∈ [0,1] constraint, it is therefore desirable to which each solution is judged.If taken except constraint
Control point gray value λ i can be obtained in boundary valuejThen k enters step S5.
Step S5 establishes the skeleton model of bone according to the parameter value of sampled point and gray value, then, enters step S6；
Step S6, to obtaining emulation strain analysis figure after the geometrical analysis such as skeleton model carries out.
In the present embodiment, body parameterized model provides not only surface information, also within control point form given
The geometry carrier of material information in body domain.Given more control point naturally also can more accurately express material information.If should
For model for analyzing, mesh refinement is also the important means for improving Finite element analysis results accuracy, therefore studies body parametrization
The refinement of model is necessary.Geometric expression precision determined by granule surface contral point, and the precision of material information expression then by
All control points codetermine, wait geometrical analysis comprising the following steps:
Step S6-1, according to the needs of practical problem, in NURBS (B-spline) unit basis that start node generates, choosing
Suitable grid subdivision strategy to be selected grid is finely divided to obtain unit grid, the subdivision of grid will not both change geometry,
Computational accuracy can also be improved, subsequently into step S6-2.
Step S6-2 introduces equal ginseng concept, the unknown quantity that will be solved, such as temperature, displacement, stress are expressed as interpolation
Basic function basic function expression-form identical with geometrical model, wherein related coefficient is generally control point or freedom degree variable, so
After enter step S6-3.
Step S6-3, referring to finite elements " weak " solution form and minimum potential energy principal by the unit of the unit grid of generation
Matrix is attached in global stiffness matrix, subsequently into step S6-4.
Step S6-4 is handled global stiffness matrix according to given boundary and constraint condition, solves linear equation
Group solves the value of other unknown quantitys by the value at the control point solved, and geometrical analysis is waited to have using the solution of " weak solution " form
First equation is limited, obtains numerical solution by way of equivalence integral, in solution procedure, most important link is element stiffness square
Battle array is assembled into global stiffness matrix, is calculated to simplify, and replaces unit integrated with the summation of Gauss integration point, according to limited
Metatheory, the linear equation of solution are to solve linear equation:
Ku=F (10)
In formula, K is the global stiffness matrix, and u is transposed matrix, and F is loading matrix, and wherein K is by described in assembling
Element stiffness matrix KeIt obtains,
The element stiffness matrix is calculated with following formula:
For the method for asking of integral function, in parameters unit theory, each unit in physical domain should be mapped to phase
Same unit, which is the tensor product of (ξ, η, ζ).If current Gauss integration point is (ξm, ηm, ζm), then we can obtain
To the new appropriate u of nodeg,vg,wg, for the range of a unit and Gauss point, ug,vg,wgIt is commonly used according in finite element method
Gauss integration point selection rule be calculated,
In above formula (11):
Rs=Ni(ug)Nj(vg)Nk(wg) (14)
In formula, S=i+j+k,
i∈[1,p+1],j∈[1,q+1],k∈[1,r+1]s∈[1,ne],ne=(p+1) (q+1) (r+1),
In addition, De is a symmetrical matrix, for homogeneous material, each place De is identical, if set
Because geometrical analysis is waited also to need the material information at Gauss integration point, i.e. Young's modulus and Poisson's ratio, and bone
Belonging to heterogeneous material, its Young's modulus and Poisson's ratio changes with position, De can be approached there are three types of method:
1. the cell level constant based on approximating unit center usually myopia De, De be usually at the center of unit estimation and
Come；
2. the De based on approximate node is estimated from cell node；
3. based on the De at approximate Gaussian integration method being estimated from the Gauss integration point of unit；
Root is needed at Gauss integration point since we solve element stiffness matrix using the method for Gauss integration
Two variables are subjected to interpolation between point according to B-spline basic function, shown in calculation formula such as formula (19),
In formula, it is Poisson's ratio that E, which is Young's modulus, ν, subsequently into step S6-4.
Step S6-4, through the above steps, geometry, material, performance information at control point can be obtained.If it is desired to knowing
The information of other points other than road control point, it is necessary to be obtained by the interpolation at control point.These three types of information in formula (11)
It is indicated respectively with G, M, R, Pi,j,k(G, M, R) indicates to assign the control point of geometry, material, performance information, and wherein performance information is logical
The geometrical analysis such as the 4th step are crossed to acquire, and geometry is acquired with material information by first three step,
Then, the model result completed to analysis is shown.
In the present embodiment, the heterogeneous body being suitable for towards a material reverse ginseng is developed using VS2008 and OpenGL
The program of numberization model construction.CT sweeps anchor picture from the age 25 years old, height 1.73m, previously without femur history of disease, outside no leg
The adult man for hurting history, using 64 row's spiral CT of machine, operating voltage 120kV, pixel 0.43mm, layer is away from 0.4mm, and totally 516 layers
The CT image of 512x512 DICOM formats.
Point cloud data is exported the array of a n × 4 dimension by software mimics by reading CT picture, and data first three columns are a little
Coordinate, the 4th column are the gray value h of the point, and reading data amount check is 116530, and data volume is huge, appropriate to simplify, and modulus type
Top be divided into example.
Fig. 2 is the volume mesh display figure of bone modeling and simulation method in the embodiment of the present invention；
Fig. 3 is the exterior view of the skeleton model of bone modeling and simulation method in the embodiment of the present invention；And
Fig. 4 is the strain analysis figure of bone modeling and simulation method in the embodiment of the present invention.
As in Figure 2-4, geometric space is carried out to simplified model and material space is rebuild, available Fig. 3 and Fig. 4
Shown in body parameterized model, represent material distribution, Fig. 2 is that volume mesh show, shared 9*17*13=1989 control point, Fig. 3
It is shown for surface.Fig. 4 is to apply point pressure in the top of model, using etc. the strain analysis figure that obtains after geometrical analysis.From
And analyzing bone bone distribution of force situation in the collision process of simulation.
The action and effect of embodiment
The bone modeling and simulation method according to involved in the present embodiment, by human skeleton carry out CT scan and
Point cloud data is extracted after being pre-processed, and then, basal plane projection is carried out according to the profile of point cloud data and obtains point cloud data
Parameter value, also, the scattered points in point cloud data are subjected to base curved surface projection and obtain the parameter value of scattered points, then, using most
Scattered points are carried out B-spline surface fitting and obtain granule surface contral point by small square law principle；It is calculated using discrete Coons body interpolation
Method obtains body parameterized model, finds the nearest node of sampled point, initial parameter value of the parameter value of the node as sampled point,
Then near initial parameter value AMSAA discrete model region, find out the B-spline body node nearest from sampled point in the zone, will
Then the parameter value of the B-spline body node as sampled point acquires surface control by principle of least square method and quadratic form optimization
The gray value of system point and interior control point, then obtains skeleton model, finally, carrying out grid subdivision according to B-spline obtains element mesh
Unknown quantity is expressed as form corresponding with the basic function of B-spline by lattice, then uses weak form and minimum potential energy principal
The element stiffness matrix of unit grid is attached in global stiffness matrix, then global stiffness matrix is handled, is obtained
The value of unknown quantity, finally to the value of unknown quantity, treated that result is shown, so, the bone modeling and simulation of the present embodiment
Method can construct more accurate skeleton model, and energy accurate expression real material distribution, be more in line with to design
The bone product of actual conditions.
In the present embodiment, due to being pre-processed to picture, so that the skeleton model extracted is more accurate.
In the present embodiment, due to using OBBs algorithm, arithmetic speed can be accelerated in this way.
Above embodiment is preferred case of the invention, the protection scope being not intended to limit the invention.
Claims (5)
1. a kind of bone modeling and simulation method, for using human skeleton as heterogeneity product according to the material of the bone
Expect that performance carries out modeling and simulation to the bone, which comprises the following steps:
Step 1, CT scan is carried out to obtain image data to the bone, then, described image data is carried out at image
Point cloud data is obtained after reason；
Step 2, basal plane projection is carried out according to the profile of the point cloud data, and corresponding with B-spline to the point cloud data
B-spline parametrization is carried out, so that the parameter value of the point cloud data is obtained, then, by the click-through at random in the point cloud data
Row base curved surface projection obtains the parameter value of the scattered points, then, the scattered points is carried out B-spline using principle of least square method
Surface fitting and obtain granule surface contral point；
Step 3, body interpolation is carried out using the granule surface contral point as known conditions, is obtained using discrete Coons body interpolation algorithm
Body parameterized model, the granule surface contral point in the body parameterized model correspond to the surface information of the bone, the body
Interior control point in parameterized model corresponds to the material information of the bone；
Step 4, multiple sampled points are set, and the body parameterized model is marked off into multiple nodes, then, find with it is described
The nearest node of sampled point, initial parameter value of the parameter value of the node as the sampled point, then described initial
The B-spline body segment nearest apart from the sampled point is found out in the region of the discrete B-spline near parameter value in this region
Point, using the B-spline body node as the parameter value of the sampled point；Also, according to the parameter value of the sampled point and pass through
Principle of least square method and quadratic form optimization acquire the gray value of the granule surface contral point and the interior control point,
Parameterizing body Model according to heterogeneous B-spline indicates are as follows:
In formula, Ni,p,Nj,q,Nk,rFor the B-spline basic function of 3D solid in three directions；H and λijkRespectively indicate the section
The grayscale information of point and the control point, and normalized has been carried out,
If the sampled point is Ts={ xs,ys,zs,hs, i.e., each sampled point meets:
According to matrix properties, above formula can also be write as following form:
The quantity of the node is ten times of the sampled point quantity；
Step 5, the skeleton model of the bone is obtained according to the parameter value of the sampled point and the gray value；And
Step 6, grid subdivision is carried out according to the B-spline and obtains unit grid, unknown quantity is expressed as the base with the B-spline
Then the corresponding form of function uses weak form and minimum potential energy principal by the element stiffness matrix of the unit grid
It is attached in global stiffness matrix, then the global stiffness matrix is handled, obtains the value of the unknown quantity, finally to institute
Stating the value of unknown quantity, treated that result is shown,
The global stiffness matrix is handled, replaces unit integrated using the summation of Gauss integration point, solves linear side
Journey:
Ku=F,
In formula, K is the global stiffness matrix, and u is transposed matrix, and F is loading matrix, and wherein K is by assembling the unit
Stiffness matrix KeIt obtains,
The element stiffness matrix is calculated with following formula:
In above formula:
Rs=Ni(ug)Nj(vg)Nk(wg),
In formula, S=i+j+k,
i∈[1,p+1],j∈[1,q+1],k∈[1,r+1]s∈[1,ne],ne=(p+1) (q+1) (r+1),
ug,vg,wgIt is calculated according to Gauss integration point selection rule common in finite element method,
Also, it sets
Young's modulus and Poisson's ratio change with position in heterogeneous material, therefore use Gauss integration method, in Gauss integration point
The Young's modulus and the Poisson's ratio are carried out interpolation according to the B-spline basic function by place between point, i.e.,
In formula, it is the Poisson's ratio that E, which is the Young's modulus, ν,.
2. bone modeling and simulation method according to claim 1, it is characterised in that:
Wherein, the step 1 includes to carry out noise filtering to described image data, sharpen feature and the pretreatment of edge enhancing.
3. bone modeling and simulation method according to claim 2, it is characterised in that:
Wherein, the step 1 is also comprising setting the threshold value of described image data, so that described in described image data
The contrast of bone increases, and handles and define the range of threshold value by the interpolated value to described image data, and use area
Domain increases principle and is split to described image data, to obtain point cloud data described in model sum of the grayscale values.
4. bone modeling and simulation method according to claim 1, it is characterised in that:
Wherein, in the step 4, the node nearest apart from the sampled point is found using OBBs algorithm.
5. bone modeling and simulation method according to claim 1, it is characterised in that:
Wherein, using the principle of least square method, the optimization problem that following formula indicates is solved:
N is set in formulaijks=Ni,p(us)Nj,q(vs)Nk,r(ws), wherein Γ (λijk) it is to meet the principle of least square method
λijk, it is rewritten into:
α=(i-1)+(j-1) × l+ (k-1) × l × m+1 is set, optimizing expression is obtained:
According to the property of B-spline basic function, obtain:
Also, last solution are as follows:
Due to λαConstraint more than or equal to 0 and less than or equal to 1, obtains the control point gray value λijk。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510689648.9A CN105224764B (en) | 2015-10-22 | 2015-10-22 | Bone modeling and simulation method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510689648.9A CN105224764B (en) | 2015-10-22 | 2015-10-22 | Bone modeling and simulation method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105224764A CN105224764A (en) | 2016-01-06 |
| CN105224764B true CN105224764B (en) | 2019-01-29 |
Family
ID=54993730
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201510689648.9A Expired - Fee Related CN105224764B (en) | 2015-10-22 | 2015-10-22 | Bone modeling and simulation method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105224764B (en) |
Families Citing this family (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105964915B (en) * | 2016-06-02 | 2018-02-09 | 林勇 | A kind of wolf towel worm nest reticulate form forming method |
| CN106384384B (en) * | 2016-09-18 | 2020-05-05 | 上海理工大学 | Shape optimization method of three-dimensional product model |
| CN107767452A (en) * | 2017-10-10 | 2018-03-06 | 上海理工大学 | The AMF general file generation methods of heterogeneous solid parameterized model |
| CN110176296A (en) * | 2019-05-11 | 2019-08-27 | 北京工业大学 | One kind being used for human femur under loading modeling method |
| CN110833472B (en) * | 2019-11-22 | 2020-12-11 | 吉林大学 | Manufacturing method of 3D printing-based individual customized knee joint bionic prosthesis |
| CN113408174B (en) * | 2021-06-28 | 2024-07-12 | 大连理工大学 | Bone model construction method, device, computer equipment and storage medium |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104039266A (en) * | 2011-11-15 | 2014-09-10 | 特里斯佩拉牙科公司 | Method and system for acquiring data from an individual for preparing a 3D model |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20110295565A1 (en) * | 2010-05-25 | 2011-12-01 | Ozen Engineering Inc. | Methods and systems of integrated simulations for patient-specific body embedded with medical implants |
-
2015
- 2015-10-22 CN CN201510689648.9A patent/CN105224764B/en not_active Expired - Fee Related
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104039266A (en) * | 2011-11-15 | 2014-09-10 | 特里斯佩拉牙科公司 | Method and system for acquiring data from an individual for preparing a 3D model |
Non-Patent Citations (2)
| Title |
|---|
| 基于等几何分析的股骨模型静力学分析;高建伟等;《计算机仿真》;20150515;第340-343页 |
| 面向材料反求的非均质体参数化模型构建;陈龙等;《机械工程学报》;20141016;第156-164页 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105224764A (en) | 2016-01-06 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105224764B (en) | Bone modeling and simulation method | |
| US8384716B2 (en) | Image processing method | |
| Hambli et al. | Multiscale methodology for bone remodelling simulation using coupled finite element and neural network computation | |
| CN112465827B (en) | Contour perception multi-organ segmentation network construction method based on class-by-class convolution operation | |
| Stoyanova et al. | An enhanced computational method for age‐at‐death estimation based on the pubic symphysis using 3 D laser scans and thin plate splines | |
| CN104392019B (en) | High-order diffusion tensor mixed sparse imaging method for brain white matter fiber tracking | |
| CN105405167A (en) | Finite element modeling method based on complete human head | |
| CN106887000A (en) | The gridding processing method and its system of medical image | |
| JPWO2006132194A1 (en) | Information processing apparatus and information processing method, image processing apparatus and image processing method, and computer program | |
| CN113781640A (en) | Three-dimensional face reconstruction model establishing method based on weak supervised learning and application thereof | |
| CN110033519A (en) | Three-dimensional modeling method, device, system and storage medium based on Implicitly function | |
| Podshivalov et al. | Multiscale FE method for analysis of bone micro-structures | |
| Schneider et al. | MedmeshCNN-Enabling meshcnn for medical surface models | |
| Coffman et al. | OOF3D: An image-based finite element solver for materials science | |
| CN112750137B (en) | Liver tumor segmentation method and system based on deep learning | |
| CN104867104A (en) | Method for obtaining anatomical structural atlas for target mouse based on XCT image non-rigid registration | |
| Gharleghi et al. | Deep learning for time averaged wall shear stress prediction in left main coronary bifurcations | |
| CN108733952B (en) | Three-dimensional characterization method for spatial variability of soil water content based on sequential simulation | |
| CN114119872A (en) | Method for analyzing 3D printing intraspinal plants based on artificial intelligence big data | |
| CN105389821A (en) | Medical image segmentation method based on combination of cloud module and image segmentation | |
| CN113764101A (en) | CNN-based breast cancer neoadjuvant chemotherapy multi-modal ultrasonic diagnosis system | |
| Warkhedkar et al. | Material-solid modeling of human body: A heterogeneous B-spline based approach | |
| CN104915989A (en) | CT image-based blood vessel three-dimensional segmentation method | |
| CN112750110A (en) | Evaluation system for evaluating lung lesion based on neural network and related products | |
| Aydin et al. | A hybrid image processing system for X-ray images of an external fixator |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20190129 Termination date: 20211022 |