JP2019205796A - Medical image processing apparatus, medical image processing method, program, and MPR image generation method - Google Patents

Abstract

To provide a medical image processing apparatus, etc. capable of clearly visualizing a specific region.SOLUTION: A medical image processing apparatus 1 includes: an image acquisition unit 21 for acquiring a plurality of tomographic images; a geometric information acquisition unit for approximating the outline of a subject region in the tomographic image by a predetermined geometric shape for each tomographic image, and acquiring information on the geometric shape; and an MPR image generation unit 24 for correcting a signal value of each pixel of the tomographic image on the basis of the information on a geometric shape corresponding to the tomographic image, and generating an MPR image on the basis of all the tomographic images whose signal values are corrected.SELECTED DRAWING: Figure 2

Description

  The present disclosure relates to a medical image processing apparatus, a medical image processing method, a program, and an MPR image generation method.

  In medical image diagnosis, there is a scene where it is desired to obtain an MPR image suitable for observation of a specific subject (human body tissue). For example, when it is desired to observe a predetermined human body tissue such as the internal organs and brain in the chest and head, the internal organs and brain are surrounded by bones, and the bone region rather hinders diagnosis. In an MPR image generated from a CT image, bones are generally clearly depicted, and internal organs and blood vessels may be hidden. Therefore, it may be necessary to devise such as cutting for visualization.

  For example, in Patent Document 1 below, extraction and labeling of regions is performed based on signal values of CT images, and whether or not each labeled rib region is appropriate is determined in the body axis direction of a separately constructed statistical rib. A technique using a statistical database having height width information is disclosed. Thereby, the bone region (rib region) can be detected and removed. However, in the method of Patent Document 1, a determination error occurs when the rib has a congenital malformation or is accompanied by a fracture. Cannot be removed.

  In addition, a 3D mask is created by extracting the bone region to be unobserved using the region expansion method (region growing method), and the bone region to be unobserved is removed by mask processing while referring to the 3D mask. A method is disclosed (see Patent Documents 2 and 3). The region expansion method is a method of designating pixels in a non-observation target region and extracting neighboring pixels one after another using that pixel as a start point (expansion start point).

JP 2009-240569 A Japanese Patent No. 4087517 Japanese Patent No. 5257958

  However, the creation of a three-dimensional mask using the region expansion method requires a user to specify a plurality of expansion start points, and the user also has to indicate the expansion termination point. There is a problem that variations occur. In addition to the region expansion method, creation of a three-dimensional mask requires user skill and experience, and the creation takes time and effort.

  The present disclosure has been made in view of the above-described problems, and an object of the present disclosure is to provide a medical image processing apparatus that can clearly visualize a specific area.

  According to an embodiment of the present disclosure, image acquisition means for acquiring a plurality of tomographic images, and for each of the tomographic images, an outline of a subject area of the tomographic image is approximated with a predetermined geometric shape, and information on the geometric shape is acquired. And a geometric information acquisition means for correcting the signal value of each pixel of the tomographic image based on information of a geometric shape corresponding to the tomographic image, and an MPR image based on all the tomographic images whose signal values are corrected. And an MPR image generation means for generating a medical image processing apparatus.

  According to an embodiment of the present disclosure, the signal value of each pixel of the tomographic image is corrected based on geometric shape information that approximates the outline of the subject area of the tomographic image, and the signal value is corrected based on the tomographic image. To generate an MPR image. Thus, by controlling the signal value based on the subject shape (geometric shape), it is possible to easily obtain an MPR image in which a specific area is clearly visualized.

  Further, the MPR image generation means further comprises correction magnification calculation means for calculating a correction magnification corresponding to each pixel of the tomographic image for each of the tomographic images based on information on a geometric shape corresponding to the tomographic image. May multiply the signal value of each pixel of the tomographic image by the correction magnification calculated by the correction magnification calculation means. Thus, the signal value can be corrected by multiplying the signal value of each pixel of the tomographic image by the correction magnification calculated based on the geometric information.

  Further, the correction magnification calculation means further includes a table creation means for creating a correction table for storing the calculated correction magnification, and the MPR generation means is configured to output the signal value of each pixel of the tomographic image with respect to the signal value. You may make it multiply the correction magnification obtained with reference to a correction table. Thereby, the signal value can be corrected with reference to the correction table storing the correction magnification.

  Further, the correction magnification calculation means may reduce the correction magnification as the distance from the geometric center of the geometric shape increases. As a result, the signal value is attenuated as the distance from the center of the subject area becomes invisible, so that the bone area and the body surface away from the center of the subject area become invisible, and an organ or the like near the center of the subject area is clearly visualized.

  The geometric shape is an ellipse, and the geometric information acquisition means acquires, as the geometric information, the center coordinates of the ellipse, the horizontal size of the ellipse, and the vertical size of the ellipse, and the correction magnification. The calculating means may calculate the correction magnification based on the center coordinates, the horizontal size, and the vertical size. Thereby, the outline of the subject region can be approximated by an ellipse, and the correction magnification can be calculated based on the ellipse parameters (center coordinates, horizontal size, vertical size). In the present disclosure, the horizontal direction is a coordinate axis that is closer to the left-right axis direction of the subject shown in the tomographic image among the two coordinate axis directions (X-axis direction and Y-axis direction) that define the coordinate system of the tomographic image. It refers to the direction. The vertical direction is the closer to the dorsoventral axis (front-rear axis) of the subject shown in the tomographic image among the two coordinate axis directions (X-axis direction and Y-axis direction) that define the coordinate system of the tomographic image. This refers to the direction of the coordinate axis. Here, the X axis and the Y axis of the tomographic image are appropriately set, for example, when photographing with a modality. In the present disclosure, the coordinate axis that is closer to the left and right axes of the subject is the X axis, and the subject's dorsoventral axis The coordinate axis whose direction is closer to the (front-rear axis) is taken as the Y-axis. That is, in the present disclosure, the horizontal direction corresponds to the X-axis direction of the tomographic image, and the vertical direction corresponds to the Y-axis direction of the tomographic image.

  Further, the correction magnification calculation means maximizes the correction magnification in the center coordinates, and the ratio of the central coordinates to the horizontal size and the vertical size is the same on each concentric ellipse inside the ellipse. The correction magnification may be decreased as the horizontal size and the vertical size of each concentric ellipse increase. Thereby, the correction magnification inside the ellipse is suitably calculated.

  Further, the correction magnification calculation means corrects the center coordinate by adding a predetermined offset to the center coordinate, and the magnification and the vertical direction with respect to the predetermined horizontal size to the horizontal size and the vertical size. The horizontal size and the vertical size are corrected by multiplying the magnification with respect to the size of the image, and the correction magnification is calculated based on the corrected center coordinates, the horizontal size, and the vertical size. It may be calculated. Thereby, the correction magnification can be calculated after correcting the ellipse parameters (center coordinates, horizontal size, vertical size).

  The image forming apparatus may further include an adjusting unit that allows a user to adjust the offset, the magnification with respect to the horizontal size, and the magnification with respect to the vertical size. Thereby, the user can adjust the correction amount of the ellipse parameter.

  Further, the MPR generation unit may correct the signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and information on a slice position of the tomographic image. Accordingly, the signal value can be corrected in consideration of the slice position information (slice parameter) of the tomographic image in addition to the geometric shape information (geometric parameter) of the tomographic image.

  Correction magnification calculating means for calculating, for each tomographic image, a correction magnification corresponding to each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and information on a slice position of the tomographic image; Further, the MPR image generation means may multiply the signal value of each pixel of the tomographic image by the correction magnification calculated by the correction magnification calculation means. Thus, the signal value can be corrected by multiplying the signal value of each pixel of the tomographic image by the correction magnification calculated based on the geometric shape information and the slice position information.

  Further, the correction magnification calculation means further includes a table creation means for creating a correction table for storing the calculated correction magnification, and the MPR generation means is configured to output the signal value of each pixel of the tomographic image with respect to the signal value. You may make it multiply the correction magnification obtained with reference to a correction table. Thereby, the signal value can be corrected with reference to the correction table storing the correction magnification.

  Further, the correction magnification calculation means may decrease the correction magnification as the distance from the geometric center of the geometric shape decreases, and decrease as the slice position is located from the center to the end. As a result, the signal value that decays away from the center of the subject region is attenuated, so that the bone region and the body surface away from the central portion of the subject region become invisible, and an organ or the like near the central portion of the subject region is clearly visualized. Further, since the signal value attenuates as the slice position moves away from the center, the bone region can be effectively removed particularly in the head CT in which the bone region (parietal bone or jaw bone) occupying the subject region increases near the end of the slice. Can do.

  The geometric shape is an ellipse, and the geometric information acquisition unit acquires, as the geometric information, an ellipse center coordinate, a horizontal size of the ellipse, and a vertical size of the ellipse, and the slice position Is the total number of slices and the slice order, and the correction magnification calculation means determines the correction magnification based on the center coordinates, the horizontal size, the vertical size, the total number of slices, and the slice order. May be calculated. Accordingly, the correction magnification can be calculated based on the ellipse parameters (center coordinates, horizontal size, vertical size) and slice parameters (total number of slices, slice rank).

  Further, the correction magnification calculation means maximizes the correction magnification in the center coordinates, and the ratio of the central coordinates to the horizontal size and the vertical size is the same on each concentric ellipse inside the ellipse. The correction magnification is decreased as the horizontal size and the vertical size of each concentric ellipse are increased, and the correction magnification is decreased as the difference between the total number of slices / 2 and the slice order is increased. May be. Thereby, the correction magnification inside the ellipse is suitably calculated in consideration of the slice direction.

  Further, the correction magnification calculation means corrects the center coordinate by adding a predetermined offset to the center coordinate, and the magnification and the vertical direction with respect to the predetermined horizontal size to the horizontal size and the vertical size. The horizontal size and the vertical size are corrected by multiplying by a magnification with respect to the size of the image, the total number of slices is corrected by multiplying by a magnification with respect to the predetermined total number of slices, and the slice order is added with a predetermined slice offset. The correction magnification may be calculated based on the corrected center coordinates, the horizontal size, the vertical size, the total number of slices, and the slice order. Thus, the correction magnification can be calculated after correcting the ellipse parameters (center coordinates, horizontal size, vertical size) and slice parameters (total number of slices, slice order).

  The image processing apparatus may further include adjustment means for allowing a user to adjust the offset, a magnification with respect to the horizontal size, a magnification with respect to the vertical size, a magnification with respect to the total number of slices, and the slice offset. Thereby, the user can adjust the correction amounts of the ellipse parameters (center coordinates, horizontal size, vertical size) and slice parameters (total number of slices, slice order).

  In addition, the approximate accuracy of the geometric shape of the tomographic image is determined based on a predetermined criterion, and information on the geometric shape of the tomographic image that does not satisfy the predetermined criterion is used as information on the geometric shape of the other tomographic image that satisfies the predetermined criterion. And geometric information correction means for correcting based on the above. Thereby, it is possible to correct the information (geometric parameters) of the geometric shape of the tomographic image with poor approximation accuracy.

  In addition, the geometric shape is an ellipse, and the geometric information acquisition unit acquires, as the geometric information, the center coordinates of the ellipse, the horizontal size of the ellipse, and the vertical size of the ellipse, and corrects the geometric information. Means for determining the horizontal size of the ellipse in the tomographic image in which the horizontal size of the ellipse matches the horizontal size of the image, and the other tomographic image in which the horizontal size of the ellipse is smaller than the horizontal size of the image; You may make it correct | amend so that it may become the same as the ratio of the size of the ellipse in the horizontal direction of the image, and the size of the vertical direction. Thereby, even if the tomographic image does not fit in the image range, the geometric parameter (ellipse parameter) can be calculated satisfactorily.

  Further, a color map acquisition unit that acquires a color map that defines a correspondence relationship between the signal value, the color value, and the opacity, and the signal value of each pixel of the tomographic image is obtained by referring to the color map. Voxel creation means for creating voxel data that is a set of voxels that convert the color value and opacity according to the color value and retain the color value and opacity, and rendering means for generating a volume rendering image based on the voxel data May be provided. Thereby, a volume rendering image can be generated based on the tomographic image.

  Further, the voxel creating means may correct the opacity according to the signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image. As a result, when the voxel data is created by converting the signal value of each pixel of the tomographic image into a color value and opacity with reference to the color map, the opacity defined in the color map (opacity curve) is converted to the tomographic image. Is corrected based on information on a geometric shape approximating the outline of the subject area. Thus, by controlling the opacity based on the subject shape (geometric shape), it is possible to easily obtain a rendered image in which a specific area is clearly visualized.

  Further, the voxel creating means may further correct the opacity according to the signal value of each pixel of the tomographic image based on information on a slice position of the tomographic image. Thereby, the opacity can be corrected by further considering the slice position information of the tomographic image.

  The voxel creating means further includes color map adjusting means for adjusting the color map so as to attenuate opacity corresponding to the signal value according to a difference value between the signal value and a predetermined threshold value. Voxel data may be created by referring to the adjusted color map. Thereby, the opacity of the color map (opacity curve) can be adjusted according to the signal value.

  A color map creating means for creating a color map that defines the correspondence between the signal value of a predetermined range, the color value, and the opacity based on the correspondence between the representative signal value, the color value, and the opacity; The color map acquisition unit may acquire the created color map. Thereby, for example, a user can set a representative signal value (for example, characteristic signal values of about 10 locations), a color value and an opacity for the signal value, and a color map can be automatically created. it can.

  According to an embodiment of the present disclosure, the control unit of the computer acquires an image acquisition step of acquiring a plurality of tomographic images, and approximates the outline of the subject area of the tomographic image with a predetermined geometric shape for each of the tomographic images. A geometric information acquisition step of acquiring information of the geometric shape, and correcting the signal value of each pixel of the tomographic image based on the information of the geometric shape corresponding to the tomographic image, and all the signal values corrected An MPR image generation step for generating an MPR image based on the tomographic image is provided.

  According to an embodiment of the present disclosure, the signal value of each pixel of the tomographic image is corrected based on geometric shape information that approximates the outline of the subject area of the tomographic image, and the signal value is corrected based on the tomographic image. To generate an MPR image. Thus, by controlling the signal value based on the subject shape (geometric shape), it is possible to easily obtain an MPR image in which a specific area is clearly visualized.

  Further, according to an embodiment of the present disclosure, the computer has an image acquisition unit that acquires a plurality of tomographic images, and for each of the tomographic images, the outline of the subject area of the tomographic image is approximated by a predetermined geometric shape, and the geometric shape Geometric information acquisition means for acquiring the information of the pixel, the signal value of each pixel of the tomographic image is corrected based on information of the geometric shape corresponding to the tomographic image, and based on all the tomographic images whose signal values are corrected A program that functions as MPR image generation means for generating an MPR image is provided. A medical image processing apparatus according to an embodiment of the present disclosure can be obtained by installing a program according to an embodiment of the present disclosure on a general-purpose computer.

  Further, according to an embodiment of the present disclosure, there is provided an MPR image generation method executed by an MPR image generation apparatus that generates an MPR image based on a tomographic image, the image acquiring step for acquiring a plurality of tomographic images, and the tomographic image For each image, a geometric information acquisition step for approximating the outline of the subject area of the tomographic image with a predetermined geometric shape and acquiring information on the geometric shape, and the signal value of each pixel of the tomographic image corresponding to the tomographic image An MPR image generation method including: an MPR image generation step of generating an MPR image based on the tomographic image corrected based on information on the geometric shape to be corrected and the signal value corrected.

  According to an embodiment of the present disclosure, the signal value of each pixel of the tomographic image is corrected based on geometric shape information that approximates the outline of the subject area of the tomographic image, and the signal value is corrected based on the tomographic image. To generate an MPR image. Thus, by controlling the signal value based on the subject shape (geometric shape), it is possible to easily obtain an MPR image in which a specific area is clearly visualized.

  According to the present disclosure, a specific area can be clearly visualized.

The figure which shows the hardware constitutions of the medical image processing apparatus 1. The figure which shows the function structure of the medical image processing apparatus 1. The figure which shows the internal structure of the correction table preparation part 23. The figure explaining a body axis section, a coronal section, and a sagittal section Flowchart showing the overall operation of the medical image processing apparatus 1 Flow chart showing correction table creation processing The figure explaining the process which acquires the ellipse parameter P (z) The figure explaining correction | amendment of the ellipse parameter P (z) The figure which shows a mode that the outline of a to-be-photographed object shape is approximated by the ellipse shape for every tomographic image Do (z). The figure which shows the distribution pattern of the correction magnification S The figure which shows the relationship between the amplitude magnification of XY direction and correction magnification The figure which shows the example of a display of MPR image The figure which shows the example of a display of MPR image The figure which shows the distribution pattern of correction magnification S (x, y, z) in 2nd Embodiment. The figure which shows the relationship between the amplitude magnification of a XYZ direction, and correction magnification. The figure which shows the function structure of the medical image processing apparatus 1 in 3rd Embodiment. The flowchart which shows the whole operation | movement of the medical image processing apparatus 1a in 3rd Embodiment. The figure which shows the example of a display of the parameter adjustment screen 40 The figure which shows the function structure of the medical image processing apparatus 1b in 4th Embodiment. The figure which shows the internal structure of the rendering part 32 The flowchart which shows the whole operation | movement of the medical image processing apparatus 1b in 4th Embodiment. The figure which shows the property of the adjusted color map Cmap Flow chart showing voxel creation processing Diagram showing the outline of the ray casting method Flow chart showing rendering process by ray casting method Flowchart showing search control mask calculation processing The figure explaining search control mask calculation processing Flow chart showing opaque voxel search process Flow chart showing ray casting process Diagram explaining 3D texture mapping method Flow chart showing color correction processing Flow chart showing rendering processing by 3D texture mapping method Flow chart showing the projection screen setting process Flow chart showing the rendering process Rendered image generated by the conventional method (front) Rendered image generated by the proposed method (front) Rendered image generated by the proposed method (back) Rendered image generated by the proposed method (side view)

  Hereinafter, embodiments of the present disclosure will be described in detail based on the drawings. In the present embodiment, a case where an MPR image is generated based on a CT image taken by an X-ray CT apparatus will be described. However, the present disclosure is based on a medical image (such as an MRI image or a PET image) other than a CT image. The present invention is also applicable when generating an MPR image.

[First Embodiment]

  FIG. 1 is a block diagram showing a hardware configuration of a medical image processing apparatus 1 (MPR image generation apparatus) in the present embodiment. As shown in FIG. 1, the medical image processing apparatus 1 includes a control unit 11, a storage unit 12, a media input / output unit 13, a communication control unit 14, an input unit 15, a display unit 16, a peripheral device I / F unit 17, and the like. And a general-purpose computer connected via the bus 18. However, the present invention is not limited to this, and various configurations can be adopted depending on the application and purpose.

The control unit 11 includes a CPU (Central Processing Unit), a GPU (Graphics Processing Unit), and a ROM (Read Only).
Memory, RAM (Random Access Memory), frame memory, etc. The CPU and GPU call and execute a program stored in the storage unit 12, ROM, recording medium, etc. to a work memory area on the RAM or frame memory, and drive and control each device connected via the bus 18. The later-described processing performed by the medical image processing apparatus 1 is realized.

  The GPU and the frame memory are mounted on the video card of the medical image processing apparatus 1. The ROM is a non-volatile memory and permanently holds a computer boot program, a program such as BIOS, data, and the like. The RAM and frame memory are volatile memories, and temporarily store programs, data, etc. loaded from the storage unit 12, ROM, recording medium, etc., and work areas used by the control unit 11 to perform various processes. Is provided.

  The storage unit 12 is an HDD (Hard Disk Drive) or the like, and stores a program executed by the control unit 11, data necessary for program execution, an OS (Operating System), and the like. As for the program, a control program corresponding to the OS and an application program for causing a computer to execute processing to be described later are stored. Each of these program codes is read by the control unit 11 as necessary, transferred to the RAM and frame memory, read to the CPU and GPU, and executed as various means.

  The media input / output unit 13 (drive device) inputs / outputs data, for example, media such as a CD drive (-ROM, -R, -RW, etc.), DVD drive (-ROM, -R, -RW, etc.) Has input / output devices. The communication control unit 14 includes a communication control device, a communication port, and the like, and is a communication interface that mediates communication between a computer and a network, and performs communication control between other computers via the network. The network may be wired or wireless.

  The input unit 15 inputs data and includes, for example, a keyboard, a pointing device such as a mouse, and an input device such as a numeric keypad. An operation instruction, an operation instruction, data input, and the like can be performed on the computer via the input unit 15. The display unit 16 includes a display device such as a liquid crystal panel, and a logic circuit or the like (video adapter or the like) for realizing a video function of the computer in cooperation with the display device. The input unit 15 and the display unit 16 may be integrated like a touch panel display.

  The peripheral device I / F (Interface) unit 17 is a port for connecting a peripheral device to the computer, and the computer transmits and receives data to and from the peripheral device via the peripheral device I / F unit 17. The peripheral device I / F unit 17 is configured by USB (Universal Serial Bus), IEEE 1394, RS-232C, or the like, and usually includes a plurality of peripheral devices I / F. The connection form with the peripheral device may be wired or wireless. The bus 18 is a path that mediates transmission / reception of control signals, data signals, and the like between the devices.

  The medical image processing apparatus 1 may be configured by a single computer, or a plurality of computers may be configured via a network. For example, when the medical image processing apparatus 1 includes a server and a client terminal, data input may be received at the client terminal, the server may perform various processes, and the client terminal may display the processing results. . In the following description, an example in which the medical image processing apparatus 1 is configured by one computer will be described as a simple configuration example.

  FIG. 2 is a diagram illustrating a functional configuration of the medical image processing apparatus 1. As illustrated in FIG. 2, the medical image processing apparatus 1 includes an image acquisition unit 21, a gradation compression unit 22, a correction table creation unit 23, an MPR generation unit 24, an MPR display unit 25, and a data output unit 26.

  The image acquisition unit 21 acquires a tomographic image group Do photographed by the CT apparatus. The tomographic image group Do includes a plurality of tomographic images (CT images) obtained by continuously photographing a subject (human body) at a predetermined slice interval along the head-to-tail axis. Each tomographic image is two-dimensional image data in DICOM format. The DICOM format is an image format generally used for medical images including a header portion and an image data portion in one file, and can store parameters and diagnostic information at the time of image capturing.

  One tomographic image is, for example, an image (CT image) of 512 × 512 pixels. Each pixel of the tomographic image is given a signal value v. In the case of a CT image, the signal value v is a CT value. In this embodiment, the signal value v (CT value) is 16-bit (−32768 ≦ v ≦ 32767) data (however, the number of bits of the signal value v is not particularly limited).

  The CT value is an X-ray attenuation coefficient of tissue expressed with water as a reference, and the type of tissue or lesion can be determined by the CT value (unit: HU (Hounsfield Unit)). The CT value is standardized by the X-ray attenuation coefficient of water and air. The CT value of water is 0 and the CT value of air is -1000. The CT value of fat is about -120 to -100, the CT value of normal tissue is about 0 to 120, and the CT value of bone is about 1000.

  The tomographic image group Do is obtained by stacking tomographic images, which are XY two-dimensional data, in the Z-axis direction, and can be expressed as XYZ three-dimensional data. For example, the tomographic image group Do is defined as XYZ three-dimensional data (Do (x, y, z)) as follows. In the present disclosure, the X axis corresponds to the left and right axis of the human body, the Y axis corresponds to the dorso-abdominal axis (front-rear axis) of the human body, and the Z axis corresponds to the head-tail axis (vertical axis) of the human body.

(Formula 1)
−32768 ≦ Do (x, y, z) ≦ 32767
0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1
Resolution: Rxy, Resolution Rz

  In (Expression 1), Sx represents the number of pixels in the horizontal direction (X-axis direction) of the tomographic image, Sy represents the number of pixels in the vertical direction (Y-axis direction) of the tomographic image, and Sz represents the number of slices (the number of tomographic images). . Rxy is the resolution of the tomographic image in the X-axis direction and Y-axis direction, and indicates the reciprocal of the pixel interval, that is, the number of pixels per unit distance. Rz is the resolution of the slice, and represents the reciprocal of the slice interval (for example, 0.5 mm or 1 mm) of the tomographic image, that is, the number of slices per unit distance.

  In addition, when representing a predetermined z-th (0 ≦ z ≦ Sz−1) tomographic image in the tomographic image group Do, it may be expressed as “tomographic image Do (z)”.

  The gradation compression unit 22 performs gradation compression on the signal value of the tomographic image group Do acquired by the image acquisition unit 21 to 8 bits as necessary (see (Expression 6) described later). A tomographic image Do (x, y, z) that has been tone-compressed to 8 bits is defined as follows.

(Formula 2)
0 ≦ Do (x, y, z) ≦ 255
0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1
Resolution: Rxy, Resolution Rz

  The correction table creation unit 23 creates a correction table S that stores the correction magnification corresponding to each pixel of the tomographic image group Do, and outputs the correction table S to the voxel creation unit 31.

As shown in FIG. 3, the correction table creation unit 23 includes a geometric information acquisition unit 23a, a geometric information correction unit 23b, and a correction magnification calculation unit 23c.
For each tomographic image Do (z), the geometric information acquisition unit 23a approximates the outline of the subject area of the tomographic image Do (z) with a predetermined geometric shape, and obtains geometric information (geometric parameter P (z), which will be described later). (Refer to (Formula 4)).

  The geometric information correction unit 23b determines the approximate accuracy of the geometric shape based on a predetermined standard, and uses the geometric parameter P (z) of the tomographic image Do (z) whose approximation accuracy is not good to another tomographic image Do with good approximation accuracy. Correction is performed based on the geometric parameter P (z) of (z) (see (Expression 5) described later).

  For each tomographic image Do (z), the correction magnification calculator 23c calculates each pixel (x, y, z) of the tomographic image Do (z) based on the geometric parameter P (z) corresponding to the tomographic image Do (z). ) Is calculated, and a three-dimensional correction table S storing the correction magnification S (x, y, z) is created (see (Expression 6) described later).

  At this time, the correction magnification S (x, y, z) is the magnification as the distance from the geometric center of the geometric shape approximating the outline of the subject area of the tomographic image Do (z) (the center of the ellipse in the embodiment described later) increases. Is calculated to be small.

  The correction magnification S (x, y, z) is a three-dimensional data array and can be handled as it is as a data table (correction table), so the correction magnification and the correction table are substantially the same. For this reason, the correction magnification and the correction table are represented by the same symbol “S”.

  The MPR generation unit 24 generates and displays an image of an arbitrary cross section based on the tomographic image group Do. For example, as shown in FIG. 4, a body axis section (axial section) parallel to the XY plane, a coronal section (coronal section) parallel to the XZ plane, and a sagittal section (sagittal section) parallel to the XY plane (Sagital cross section) MPR image is generated. Further, an MPR image of an oblique section (oblique section) obtained by obliquely cutting the XYZ space may be generated.

  In particular, in the present disclosure, the MPR generation unit 24 performs the correction magnification S (x) corresponding to each pixel (x, y, z) with respect to the signal value v of each pixel (x, y, z) of the tomographic image group Do. , Y, z) is multiplied to correct the signal value v, and an MPR image is generated based on the corrected signal value v.

  As described above, the correction magnification S (x, y, z) decreases as the distance from the geometric center (center of gravity) of the geometric shape approximating the outline of the subject region decreases. The signal value is attenuated. As a result, the bone region (sternum and skull), skin, external bed, fixing jig, etc. become invisible, and an MPR image in which an organ near the center of the subject region is clearly visualized can be obtained.

  The MPR display unit 25 displays the MPR image generated by the MPR generation unit 24 on the display unit 16. Further, in order to compare the results before and after signal value correction, MPR images generated based on the tomographic image group Do (original image) before correction may be displayed together.

  The data output unit 26 outputs various data and stores the data in the storage unit 12 by receiving a data saving operation from the user. For example, the data output unit 26 outputs and stores the MPR image generated and displayed by the MPR generation unit 24 and the MPR display unit 25.

Next, the operation of the medical image processing apparatus 1 will be described with reference to FIGS.
FIG. 5 is a flowchart showing the overall operation of the medical image processing apparatus 1.

  First, the control unit 11 (image acquisition unit 21) of the medical image processing apparatus 1 acquires a plurality of tomographic images (tomographic image groups) Do in the DICOM format, which are MPR image generation targets (step S1 in FIG. 5). Although the acquisition method of the tomographic image group Do is arbitrary, for example, when the tomographic image group Do is stored in the storage unit 12 in advance, the control unit 11 reads and acquires the tomographic image group Do from the storage unit 12. Can do.

  When the tomographic image group Do is stored in a storage medium such as a USB, the control unit 11 can read and acquire the tomographic image group Do from the storage medium connected to the peripheral device I / F unit 17. . Alternatively, when the tomographic image group Do is stored in the server, the control unit 11 can download and acquire the tomographic image group Do from the server via the communication control unit 14.

  Subsequently, the control unit 11 (gradation compression unit 22) performs gradation compression processing on the tomographic image group Do as necessary (step S2 in FIG. 5). Specifically, the control unit 11 converts the 16-bit tomographic image group Do (see (Equation 1)) into an 8-bit tomographic image group Do (see (Equation 2)) in the following procedure.

(Formula 3)
(1) The minimum value Dmin and the maximum value Dmax of all pixels in the Sz / 2-th intermediate slice Do (z / 2) of the tomographic image group Do are calculated.
(2) The lower limit value Lmin = (Dmax−Dmin) · γ + Dmin and the upper limit value Lmax = (Dmax−Dmin) · (1−γ) + Dmin are set. Here, γ is the contrast adjustment width of the gradation-compressed image, and the contrast increases as it approaches 0 (however, the luminance decreases). Normally, γ = 0.1 is set.
(3) The tomographic image group Do is converted into an 8-bit tomographic image group Do by the following calculation formula (compressed to 256 levels within the range from the lower limit value Lmin to the upper limit value Lmax).
Do (x, y, z)
= (Do (x, y, z) -Lmin) .255 / (Lmax-Lmin)
However, when Do (x, y, z)> 255, Do (x, y, z) = 255, and when Do (x, y, z) <0, Do (x, y, z) = 0. Saturate.

Subsequently, the control unit 11 (correction table creating unit 23) creates a correction table S storing the correction magnification corresponding to each pixel of the tomographic image group Do (step S3 in FIG. 5).
Here, the process of creating the correction table S in step S3 will be described with reference to FIGS.

FIG. 6 is a flowchart showing a process for creating the correction table S.
First, the control unit 11 (geometric information acquisition unit 23a) approximates the outline of the subject shape of the tomographic image Do (z) with a predetermined geometric shape for each tomographic image Do (z) (0 ≦ z ≦ Sz−1). Then, the geometric parameter P (z) is acquired (step S21 in FIG. 6). In the present embodiment, the outline of the subject shape is approximated by an ellipse, and an ellipse parameter (ellipse parameter P (z)) is acquired as the geometric parameter P (z).

FIG. 7 is a diagram for explaining processing for acquiring the elliptic parameter P (z) of the tomographic image Do (z).
As illustrated in FIG. 7A, the control unit 11 fills the tomographic image Do (z) with pixels that have a signal value equal to or greater than a predetermined threshold value “white” and pixels that are less than the threshold value “black”. Binarize and extract the subject area. In the example of FIG. 7A, the “white” area of the binarized tomographic image Do (z) is the subject area. For example, when extracting a subject area other than air, “−700” is set when the signal value v is 16 bits, and “10” is set when the signal value v is 8 bits.

  Subsequently, as shown in FIG. 7B, the control unit 11 specifies the minimum rectangular area R (z) including all the extracted subject areas (white areas), and as shown in FIG. An ellipse E (z) inscribed in the rectangular region R (z) is obtained as an ellipse that approximates the outline of the subject shape. Then, as shown in FIG. 7C, the control unit 11 uses the center coordinates (Cx (z), Cy (z)) of the ellipse E (z) as the ellipse parameter P (z) of the ellipse E (z). The horizontal size W (z) and the vertical size H (z) are calculated as follows. As shown in FIG. 7C, the minimum coordinates of the rectangular region R (z) are (Xmin, Ymin) and the maximum coordinates are (Xmax, Ymax).

(Formula 4)
Center coordinates Cx (z) = (Xmax + Xmin) / 2
Cy (z) = (Ymax + Ymin) / 2
Horizontal size W (z) = Xmax−Xmin + 1
Vertical size H (z) = Ymax−Ymin + 1

  The control unit 11 obtains an ellipse E (z) that approximates the outline of the subject region for all tomographic images Do (z) (0 ≦ z ≦ Sz−1), and calculates an ellipse parameter P (z).

  Subsequently, the control unit 11 (geometric information correction unit 23b) corrects the geometric parameter P (z) (elliptical parameter P (z)) for the tomographic image Do (z) whose approximation accuracy of the geometric shape (ellipse) is not good. (Step S22 in FIG. 6).

  FIG. 8 is a diagram for explaining correction of the elliptic parameter P (z) of the tomographic image Do (z). A tomographic image Do (z) in FIG. 8A is an example of a tomographic image in which the elliptical approximation cannot be satisfactorily performed. In the tomographic image in FIG. 8A, both ends of the subject area (width direction) are cut off, and the subject area does not fall within the image range (such an image may be obtained depending on conditions at the time of photographing). In such a case, the subject area cannot be accurately approximated to an ellipse.

  That is, as shown in FIG. 8B, the horizontal size W (z) of the ellipse E (z) obtained from the tomographic image is the same as the horizontal size Sx of the tomographic image (W (z) = Sx) does not accurately reflect the actual width of the subject area (width greater than Sx). For this reason, in the present embodiment, the lateral size W (z) is corrected for the tomographic image Do (z) in which both ends of the subject region (width direction) are cut.

  Specifically, the control unit 11 first displays one or more tomographic images Do (z) in which the horizontal size W (z) of the ellipse acquired in step S21 is W (z) = Sx, Is specified as a tomographic image in which both ends of the image are cut off (a tomographic image in which the approximation accuracy of the geometric shape is not good).

  Next, for each specified tomographic image Do (z), the control unit 11 has the slice position closest to the tomographic image Do (z) and the horizontal size W (z) of the ellipse acquired in step S21. ) Specifies a tomographic image satisfying W (z) <Sx (denoted as “tomographic image Do (z ′)”).

  In the tomographic image Do (z ′), the horizontal size W (z ′) of the ellipse E (z ′) is smaller than the image size Sx, and therefore the tomographic image (geometrical) in which the subject region (width direction) is within the image range. A tomographic image with good shape approximation accuracy). Then, the control unit 11 calculates the horizontal size W (z) of the ellipse E (z) obtained from the tomographic image Do (z) with both ends of the subject region cut from the tomographic image Do (z ′). Correction is performed by the following expression so that the ratio of the horizontal size W (z ′) and the vertical size H (z ′) of E (z ′) is the same.

(Formula 5)
W (z) = H (z) · W (z ′) / H (z ′)

  Thus, as shown in FIG. 8C, the horizontal size W (z) is corrected to the horizontal size W (z) of the ellipse E (z) that approximates the outline of the true subject area. Even when the subject region (width direction) does not fit within the image range, the ellipse parameter P (z) can be calculated satisfactorily.

  As described above, by the processes in steps S21 and S22 in FIG. 6, the outline of the subject shape is approximated by an elliptical shape for each tomographic image Do (z) as shown in FIG. 9, and all tomographic images Do (z) (0 An elliptic parameter P (z) (0 ≦ z ≦ Sz−1) is obtained for ≦ z ≦ Sz−1).

Then, the control unit 11 (correction magnification calculating unit 23c), for each tomographic image Do (z), based on the geometric parameter P (z) (elliptical parameter P (z)) corresponding to the tomographic image Do (z). A correction magnification S (x, y, z) corresponding to each pixel (x, y, z) of the tomographic image Do (z) is calculated, and a correction table S is created (step S23 in FIG. 6).
For example, the correction magnification S (x, y, z) of each pixel (x, y, z) of the zth tomographic image Do (z) is calculated as follows.

(Formula 6)
S (x, y, z) = 1−k · r (x, y, z) amp
However, 0 ≦ S (x, y, z) ≦ 1
r (x, y, z) = {(x−Cx (z) −xoffset) / (a · W (z) / 2)} ² + {(y−Cy (z) −yoffset) / (b · H (Z) / 2)} ²

In (Expression 6), k is an attenuation coefficient, and the initial value is 1.0. When k = 1.0, the bones and the bed become invisible. a is a correction coefficient for the horizontal size W (z), and b is a correction coefficient for the vertical size H (z). If not corrected, a = b = 1.0 is set.
xoffset and yoffset are offset values in the XY plane, and the initial value is xoffset = yoffset = 0 (unit: voxel).

  FIG. 10 is a diagram showing a distribution pattern (concentric elliptic distribution) of the correction magnification S (x, y, z) calculated by (Equation 6). As shown in the figure, the correction magnification S at the center of the ellipse (Cx (z), Cy (z)) is 1 (maximum), and the correction magnification S decreases as the distance from the center increases. Specifically, the correction magnification on each concentric ellipse in which the ratio of the center of the ellipse (Cx (z), Cy (z)) to the horizontal size W (z) and the vertical size H (z) is the same. S decreases as the diameter (lateral size and vertical size) of each concentric ellipse increases (decreases in inverse proportion to the square of the diameter). The correction magnification S outside the ellipse is uniformly 0.

  In addition, amp in (Expression 6) is an amplitude magnification on the XY plane, and an initial value is amp = 1.0. As shown in FIG. 11, when it is desired to flatten the correction magnification at the center of the subject and widen the saturation region, amp> 1.0 is set (however, shading is weakened).

  The control unit 11 corrects the magnification S (x, y, z) (0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1) for all tomographic images Do (z). And a three-dimensional correction table S for storing the calculated correction magnification S (x, y, z) is created. As described above, the correction magnification S (x, y, z) is a three-dimensional data array and can be handled as it is as a data table (correction table), so the correction magnification and the correction table are substantially the same.

The processing for creating the correction table S has been described above with reference to FIGS.
Returning to the flowchart of FIG.

  Subsequently, the control unit 11 (MPR generation unit 24, MPR display unit 25) generates an MPR image and displays it on the display unit 16 (step S4 in FIG. 5). Specifically, the MPR image is generated and displayed after correcting the signal value v (= Do (x, y, z)) of the tomographic image group Do as follows.

  When the signal value of the tomographic image group Do is 16 bits (see (Equation 1)), the control unit 11 performs the following procedure to obtain the signal value v (= Do (x, y, z)) of the tomographic image group Do. Correct.

(Formula 7)
(1) In the tomographic image group Do, the minimum value Dmin and the maximum value Dmax of all pixels in the Sz / 2-th intermediate slice Do (z / 2) are calculated.
(2) The lower limit value Lmin = (Dmax−Dmin) · γ + Dmin is set. Here, γ is the contrast adjustment width of the gradation-compressed image. The closer the value is to 0, the greater the contrast but the lower the luminance. Normally, γ = 0.1 is set.
(3) The signal value v (= Do (x, y, z)) of the tomographic image group Do is corrected as follows.
If Do (x, y, z) ≧ Lmin:
Do (x, y, z) = (Do (x, y, z) −Lmin) · S (x, y, z) + Lmin
When Do (x, y, z) ≦ Lmin:
Do (x, y, z) = Lmin

  On the other hand, when the signal value of the tomographic image group Do is 8 bits (see (Equation 2)), the control unit 11 corrects the signal value v (Do (x, y, z)) as follows.

(Formula 8)
Do (x, y, z) = Do (x, y, z) · S (x, y, z)

  Do (x, y, z) of (Expression 7) = (Do (x, y, z) −Lmin) · S (x, y, z) + Lmin, and Do (x, y, z of (Expression 8)) ) = Do (x, y, z) · S (x, y, z) is a process according to the gist of the present disclosure, whereby the correction magnification S with respect to the signal value of each pixel of the tomographic image group Do. Is multiplied to correct the signal value. In particular, the correction magnification S (x, y, z) decreases as the distance from the center of the ellipse approximating the outline of the subject region decreases, so that the signal value decreases as the pixel is farther from the center of the subject region.

  12 and 13 are diagrams illustrating display examples of MPR images generated based on the head CT image. 12 and 13 show the MPR image generated based on the tomographic image group Do (original image) before correction and the MPR image generated based on the tomographic image Do after correction. It is the example comparatively displayed about each cross section of a cross section and a sagittal cross section. As shown in the figure, by correcting the signal value v with the correction magnification S, the signal value in the skull region is significantly attenuated.

  The control part 11 (data output part 26) can output various data and can preserve | save in the memory | storage part 12 by receiving the preservation | save operation of data from a user (step S5 of FIG. 5). For example, the control unit 11 outputs and stores the MPR image (three-dimensional image) generated and displayed in step S4.

  Heretofore, the first embodiment of the present disclosure has been described. According to the first embodiment, the signal value v of each pixel (x, y, z) of the tomographic image group Do is based on ellipse information (ellipse parameter P (z)) approximating the outline of the subject area. The signal value is corrected by multiplying the calculated correction magnification S (x, y, z) (see (Expression 7) and (Expression 8)). Then, an MPR image is generated and displayed based on the tomographic image group Do whose signal values are corrected. Thereby, an MPR image in which a specific area is clearly visualized can be obtained without creating a three-dimensional mask as in the prior art. In this disclosure, it is necessary to create a correction table related to the correction magnification S (x, y, z) in place of the conventional three-dimensional mask. However, these methods do not take a method of creating a correction table in advance as described above. In addition, it is possible to use a method of calculating for each pixel of the tomographic image group Do each time according to (Equation 6), and it is not necessary to ensure on the memory as in the case of a three-dimensional mask.

  The correction magnification S decreases as the distance from the center of the ellipse approximating the outline of the subject region decreases, so that the signal value decreases as the distance from the center of the subject region increases, and the bone region (sternum or skull) or Skin, an external bed, a fixing jig, etc. become invisible, and an organ near the center of the subject area is clearly visualized.

  The correction magnification S is analytically calculated based on an ellipse parameter P (z) (center coordinate, horizontal size, vertical size) that can express the outline shape of the subject region with a small number of variables. Therefore, the signal value can be controlled without complicated processing.

  In addition, the ellipse parameter P (z) of the tomographic image Do (z) whose ellipse approximation accuracy is not good is corrected based on the ellipse parameter P (z) of another tomographic image Do (z) with good approximation accuracy. Thus, the ellipse parameter P (z) can be calculated satisfactorily even when the ellipse approximation cannot be satisfactorily performed due to shooting conditions or the like.

[Second Embodiment]
Since the head is nearly spherical, the area ratio of the bone region occupying the subject region is greatly changed in the slice direction (Z-axis direction) as compared to the chest. In particular, the area ratio of the bone area (parietal bone and jawbone) occupying the subject area increases near the end of the slice. It may be difficult to visualize. For example, since the central region of the tomographic image is also a bone region at the top of the head, the bone region cannot be made invisible with a concentric ellipse pattern. Therefore, in the second embodiment, the correction magnification S is calculated further considering the slice position of the tomographic image, thereby improving the accuracy of the MPR image particularly in the head CT.

In the second embodiment, the process in step S23 of FIG. 6 in the first embodiment (process for creating the correction table S) is changed as follows. That is, in step S23 of FIG. 6, the control unit 11 (correction magnification calculation unit 23c) for each tomographic image Do (z), the geometric parameter P (z) (elliptical parameter P ( z)) and slice position information of the tomographic image Do (z) (hereinafter referred to as “slice parameter Ps (z)”), corresponding to each pixel (x, y) of the tomographic image Do (z). A correction magnification S (x, y, z) is calculated, and a three-dimensional correction table S in which the calculated correction magnification S (x, y, z) is stored is created.
For example, the correction magnification S (x, y, z) of each pixel (x, y) of the z-th tomographic image Do (z) is calculated as follows.

  First, similarly to (Expression 6) in step S23 of FIG. 6, the control unit 11 corrects the correction magnification based on the geometric parameter P (z) (ellipse parameter P (z)) corresponding to the tomographic image Do (z). S (x, y, z) is calculated as follows.

(Formula 9)
S (x, y, z) = 1−k · r (x, y, z) · amp
However, 0 ≦ S (x, y, z) ≦ 1
r (x, y, z) = {(x−Cx (z) −xoffset) / (a · W (z) / 2)} ² + {(y−Cy (z) −yoffset) / (b · H (Z) / 2)} ^2s

  Furthermore, in the second embodiment, the control unit 11 calculates the correction magnification Sz (z) in the Z direction based on the slice parameter Ps (z) as follows, and S (x, y, z) Multiply by. Here, the slice parameter Ps (z) refers to the total number of slices Sz and the slice order z.

(Formula 10)
Sz (z) = 1.0−kz · {(z−Sz / 2−zoffset) / Sz · c / 2}} ² · ampz
S (x, y, z) = S (x, y, z) · Sz (z)
However, 0 ≦ S (x, y, z), Sz (z) ≦ 1

In (Expression 9) and (Expression 10), k is an XY plane, and kz is an attenuation coefficient in the Z direction. The initial values are k = 1.0 and kz = 0.0. In this case, the bone and the bed are invisible. In the case of the head CT, kz> 0.0 is set.
a and b are correction coefficients of the horizontal size W (z) and the vertical size H (z), respectively, and c is a correction coefficient (slice total magnification) of the Z-direction width Sz / 2. If not corrected, a = b =
Set c = 1.0.
xoffset and yoffset are XY planes, zoffset is an offset value in the Z direction, and an initial value is xoffset = yoffset = 0 (unit: voxel).

  According to (Equation 10), as the slice position of the tomographic image Do (z) is located from the center to the end (as the difference between the total number of slices Sz / 2 and the slice order z increases), the correction magnification S (x, y, z) is uniformly attenuated, and S (x, y, z) takes a value of less than 1 even at the center of the ellipse. As a result, the signal value of the bone region (parietal bone or jaw bone) near the end of the slice can be effectively attenuated.

FIG. 14 is a diagram showing a distribution pattern of the correction magnification S (x, y, z).
The upper diagram of FIG. 14 shows the distribution pattern (concentric ellipse distribution) of the correction magnification S (x, y, z) based only on the ellipse parameter P (z) calculated by (Equation 9). Indicates a distribution pattern (concentric elliptic distribution) of the correction magnification S (x, y, z) based on the ellipse parameter P (z) and the slice parameter Ps (z) calculated by (Equation 10). As shown in FIG. 14, by multiplying by the Z-direction correction magnification Sz (z), the correction magnification S (x, y, z) is uniformly attenuated according to the slice position. That is, the distribution pattern of the correction magnification S is a two-dimensional concentric ellipsoid distribution in the first embodiment, but is expanded to a three-dimensional concentric ellipsoid distribution in the second embodiment.

  In addition, amp in (Expression 9) is the XY plane, ampz in (Expression 10) is the amplitude magnification in the Z direction, and the initial value is amp = ampz = 1.0. As shown in FIG. 15, when the correction magnification at the center of the subject is flattened and the saturation region is to be widened, amp and ampz> 1.0 are set (however, shading is weakened).

  The control unit 11 calculates the correction magnification S (x, y, z) (0 ≦ x) obtained by multiplying the correction magnification Sz (z) in the Z direction by (Equation 9) and (Equation 10) for all tomographic images Do (z). ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1) is calculated, and a three-dimensional correction table S for storing the calculated correction magnification S (x, y, z) is created. .

  As described above, according to the second embodiment, the correction magnification S is calculated in consideration of the slice position of the tomographic image, so that more flexible signal value control can be performed. This is particularly effective when generating an MPR image based on a head CT image. Since the head is close to a sphere, the area ratio of the bone region occupying the subject region changes greatly in the slice direction (Z-axis direction) compared to the chest, and especially the parietal bone near the end of the slice where the area ratio of the bone region increases In some cases, it is difficult to accurately visualize the parietal region and the brain stem region of the brain by attenuating signal values of the region of the jaw and the jawbone. According to the second embodiment, as the slice position of the tomographic image is located from the center to the end, the correction magnification S corresponding to all the pixels of the tomographic image is uniformly reduced. Even in the head where the occupying bone area (the parietal bone or the jaw bone) increases, the signal value of the bone area can be effectively attenuated.

[Third Embodiment]
In the third embodiment, the user can adjust various parameters when creating the correction table S on the GUI.

  FIG. 16 is a diagram illustrating a functional configuration of the medical image processing apparatus 1a according to the third embodiment. As shown in FIG. 16, the third embodiment further includes a parameter adjustment unit 27 in addition to the functional configuration of the medical image processing apparatus 1 in the first embodiment (see FIG. 2).

  The parameter adjustment unit 27 displays the parameter adjustment screen 40 and allows the user to adjust various parameters when the correction magnification S is calculated. Specifically, when the correction magnification S is calculated based on the ellipse parameter P (z) as shown in (Equation 6) (in the case of the first embodiment), the size of the ellipse in the horizontal direction is calculated. Correction coefficient a, vertical size correction coefficient b, X-direction offset xoffset, Y-direction offset yoffset, and XY plane attenuation coefficient k.

  On the other hand, when the correction magnification S is calculated based on the ellipse parameter P (z) and the slice parameter Ps (z) as shown in (Equation 9) and (Equation 10) (in the case of the second embodiment). In addition to the above parameters, the slice total magnification c, the Z-direction offset zoffset, and the Z-direction attenuation coefficient kz can be further adjusted.

  FIG. 17 is a flowchart showing the overall operation of the medical image processing apparatus 1a according to the third embodiment. The process of generating and displaying the MPR image (Steps S31 to S34) is the same as Steps S1 to S4 of the first embodiment in FIG. In the third embodiment, after the MPR image is generated and displayed on the display unit 16 in step S34, when the user determines that the parameter needs to be adjusted (step S35; “Yes”), the parameter adjustment screen 40 displays the parameter. Is adjusted by the user (step S36).

  FIG. 18 shows an example of the parameter adjustment screen 40. As shown in the figure, in the parameter adjustment screen 40, the XY plane attenuation coefficient k (XY direction attenuation [%] in the figure), the Z direction attenuation coefficient kz (Z direction attenuation [%] in the figure), and the horizontal size are shown. Correction coefficient a (elliptical magnification X in the figure), vertical size correction coefficient b (elliptical magnification Y in the figure), total slice number magnification c (elliptical magnification Z in the figure), offset offset xoffset in the X direction (center of the figure) Offset X), Y direction offset yoffset (center offset Y in the figure), and Z direction offset zoffset (center offset Z in the figure) can be adjusted.

  When the parameter adjustment is completed (when the “OK” button on the parameter adjustment screen 40 is pressed), the process returns to step S33, and the control unit 11 executes the processes of steps S33 to S34 again. That is, the control unit 11 creates the correction table S again based on the adjusted parameters (step S33), and generates and displays an MPR image (step S34). The processing of step S33 to step S34 is repeatedly executed until the user determines that the parameter adjustment is sufficient (step S35; “No”) while checking the MPR image generated and displayed in step S34. When a good MPR image result is obtained and the user determines that the parameter adjustment is sufficient (step S35; “No”), the control unit 11 outputs data such as an MPR image as necessary (step S37). The process is terminated. The data output process in step S37 is the same as step S5 in the first embodiment in FIG.

  As described above, according to the third embodiment, the user can adjust various parameters when creating the correction table S on the GUI. In particular, when a shooting jig (bed, fixture, etc.) is reflected in the subject area, the diameter of the ellipse that approximates the outline of the subject area is calculated to be larger than the actual subject width (head width or chest width) The center position of the ellipse may shift, and the approximation accuracy of the ellipse may deteriorate. According to the third embodiment, the ellipse diameter, center position, and the like can be adjusted as appropriate on the GUI, so that a good MPR image can be obtained even in the above case.

  In the present embodiment, an example has been described in which parameters can be adjusted after displaying the MPR image so that the parameters can be adjusted while viewing the MPR image. However, the correction table S is created not only after the MPR image is displayed. It can be adjusted at any stage prior to. Further, the correction magnification S is calculated based on (Equation 9) and (Equation 10) at the stage of regeneration of the MPR image in step S34 without re-creating the correction table S based on the adjusted parameter in step S33. A method of calculating each pixel of the group Do every time can also be used.

[Fourth Embodiment]
The fourth embodiment is further characterized in that a volume rendering image is generated and displayed.

  FIG. 19 is a diagram illustrating a functional configuration of the medical image processing apparatus 1b according to the fourth embodiment. As shown in FIG. 19, in the fourth embodiment, in addition to the functional configuration of the medical image processing apparatus 1a in the third embodiment (see FIG. 16), an area designating unit 28, a color map acquisition unit 29, a color The map adjustment unit 30, the voxel creation unit 31, and the rendering unit 32 are further provided.

The region specifying unit 28 receives a region of interest (ROI) from the user.
(Interest) is accepted. For example, when the region of interest ROI is designated by a rectangular parallelepiped, the X direction ROI (X direction open section), the Y direction ROI (Y direction open section), and the Z direction ROI (Z direction open section) are as follows. Is specified.

(Formula 11)
X direction ROI = (Xs, Xe)
Y direction ROI = (Ys, Ye)
Z direction ROI = (Zs, Ze)
When the region of interest is not designated, Xs = 0, Xe = Sx-1, Ys = 0, Ye = Sy-1, Zs = 0, Ze = Sz-1.

  The color map acquisition unit 29 acquires a color map Cmap to be applied to the tomographic image group Do. The color map acquisition unit 29 outputs the acquired color map Cmap to the color map adjustment unit 30 when further adjusting, and outputs the color map Cmap to the voxel creation unit 31 when the color map Cmap is not adjusted.

  The color map Cmap defines the correspondence between the signal value v, the color value (specifically, RGB value), and the opacity (α value), and the signal value v is converted into a 24-bit color value (RGB value). And a function (substantially a two-dimensional data table) for converting to 8-bit opacity (α value). For example, a color map Cmap applied to a 16-bit tomographic image group Do (see (Expression 1)) is defined as follows.

(Formula 12)
0 ≦ Cmap (v, n) ≦ 255
-32768 ≦ v ≦ 32767
n = 0 (R), 1 (G), 2 (B), 3 (α)

  The color map Cmap (v, n) applied to the 8-bit tomographic image group Do (see (Expression 2)) is defined as follows.

(Formula 13)
0 ≦ Cmap (v, n) ≦ 255
0 ≦ v ≦ 255
n = 0 (R), 1 (G), 2 (B), 3 (α)

  Among the color maps Cmap (v, n) (0 ≦ n ≦ 3) shown in (Expression 12) and (Expression 13), in particular, the color map Cmap (v, n) (0 ≦ n ≦ 2) has a signal value v. Corresponds to a function for converting the color value into a color value (RGB value). In the present disclosure, a color map Cmap (v, n) (0 ≦ n ≦ 2) that converts the signal value v into a color value is referred to as a “color palette”.

  Of the color maps Cmap (v, n) (0 ≦ n ≦ 3) shown in (Equation 12) and (Equation 13), in particular, the color map Cmap (v, 3) converts the signal value v into opacity. It is generally called an “opacity curve”.

  In the color map Cmap (v, 3) (opacity curve), for example, the opacity of a normal tissue (16-bit signal value v = 0 to about 120) is set to 255 (opacity = 255 is opaque (all light is reflected) ), And the opacity of the bone region (16-bit signal value v = around 1000) is set to an intermediate value of 1 to 254 (opacity = 1 to 254 is translucent (part of incident light) Is reflected and others are transmitted), and the opacity is set to 0 for other tissues and the like (opacity = 0 is transparent (all incident light is transmitted)) Is done.

  The color map adjusting unit 30 signals Cmap (v, 3) (also referred to as an opacity curve) that defines the correspondence between the signal value v and the opacity in the color map Cmap (v, n) as necessary. Adjustment is performed according to the value (see (Expression 14) and (Expression 15) described later), and the adjusted color map Cmap is output to the voxel creation unit 31.

  The voxel creation unit 31 includes the tomographic image group Do, the color map Cmap acquired by the color map acquisition unit 29 (or the color map Cmap adjusted by the color map adjustment unit 30), and the correction table created by the correction table creation unit 23. Based on S and the region of interest ROI specified by the region specifying unit 28, a voxel structure V (voxel data) is created and output to the rendering unit 32.

  Specifically, the voxel creation unit 31 refers to the color map Cmap, and converts the signal value v of each pixel (x, y, z) of the tomographic image group Do into a color value and an undefined value corresponding to the signal value v. A voxel structure V (voxel data) that is a set of voxels that convert to transparency and retain color values and opacity is created.

  In particular, in the present embodiment, the voxel creation unit 31 converts the signal value v of each pixel (x, y, z) of the tomographic image group Do into a color value and opacity by referring to the color map Cmap to generate a voxel structure. When creating the body V (voxel data), the opacity defined in the color map Cmap (v, 3) (opacity curve) is not used as it is, but the correction magnification corresponding to each pixel (x, y, z). Opacity is corrected by multiplying S (x, y, z) (see (Equation 17) described later).

  As described above, the correction magnification S (x, y, z) decreases as the distance from the geometric center (center of gravity) of the geometric shape approximating the outline of the subject region decreases. Opacity is attenuated. As a result, the bone region (the sternum and skull), the skin, the external bed, the fixing jig, and the like are made transparent, and a rendering image in which an organ near the center of the subject region is clearly visualized can be obtained.

  The rendering unit 32 executes volume rendering processing based on the voxel structure V (voxel data) created by the voxel creation unit 31, and generates and displays a volume rendering image (three-dimensional image). The rendering unit 32 generates, as a volume rendering image, a virtual rendering image for observing the bronchi, the large intestine, and the like by freely moving the entire rendering image for observing the subject from a distal viewpoint and the viewpoint to the subject.

As shown in FIG. 20, the rendering unit 32 includes a first rendering unit 32a and a second rendering unit 32b.
The first rendering unit 32a is a rendering unit that executes rendering processing by the CPU, and the second rendering unit 32b performs rendering processing by the GPU using the OpenGL (Open Graphics Library) that is an open standard of computer graphics API. A rendering unit to be executed. The rendering process is executed by the first rendering unit 32a or the second rendering unit 32b.

  FIG. 21 is a flowchart showing the overall operation of the medical image processing apparatus 1b according to the fourth embodiment. Since the processing of step S41 to step S46 is the same as step S31 to step S36 of the third embodiment in FIG. 17, it will be described from step S47 onward.

  The control unit 11 (region specifying unit 28) receives the specification of the region of interest ROI of the subject from the user as necessary (step S47 in FIG. 21). For example, the region of interest ROI is specified by a rectangular parallelepiped or the like, as described above (Formula 11).

  Subsequently, the control unit 11 (color map acquisition unit 29) acquires a color map Cmap to be applied to the tomographic image group Do (step S48 in FIG. 21). Here, when the signal value of the tomographic image group Do is 16 bits (see (Equation 1)), a 16-bit color map Cmap (see (Equation 12)) is acquired, and the signal value of the tomographic image group Do is When compressed to 8 bits (see (Expression 2)), a color map Cmap (see (Expression 13)) corresponding to 8 bits is acquired.

  The acquisition method of the color map Cmap is not particularly limited. For example, the control unit 11 allows the user to select a color map Cmap to be applied to the tomographic image group Do from a plurality of color maps Cmap prepared in advance. Can be acquired. Further, the control unit 11 may read the color map Cmap used at the time of previous rendering from the storage unit 12 and acquire it as a color map Cmap to be applied to the tomographic image group Do.

  Further, the control unit 11 may acquire a color map Cmap created by the user on a color map creation screen (not shown). In this case, the control unit 11 displays, on the color map creation screen, representative signal values v (for example, characteristic signal values at about 10 locations), color values (RGB values) for the signal values v, and opacity (α Value) is set by the user, and a color map Cmap for converting the signal value v in a predetermined range into the color value and the opacity based on the correspondence between the representative signal value v set by the user and the color value and the opacity. Generate automatically.

  Subsequently, when the color map Cmap acquired in step S48 is further adjusted (step S49; “Yes”), the control unit 11 proceeds to step S50 and adjusts the color map Cmap. Specifically, the control unit 11 (color map adjustment unit 30), as described below, out of the color map Cmap (v, n), a color map Cmap (defining the correspondence between the signal value v and the opacity. v, 3) Adjust the opacity curve.

  When the signal value of the tomographic image group Do is 16 bits (see (Expression 1)), the control unit 11 performs a 16-bit color map Cmap (v, 3) (see (Expression 12)) by the following procedure. Adjust.

(Formula 14)
(1) In the tomographic image group Do, the minimum value Dmin and the maximum value Dmax of all pixels in the Sz / 2-th intermediate slice Do (z / 2) are calculated.
(2) The lower limit value Lmin = (Dmax−Dmin) · γ + Dmin and the upper limit value Lmax = (Dmax−Dmin) · (1−γ) + Dmin are set. Here, γ is the contrast adjustment width of the gradation-compressed image, and the contrast increases as it approaches 0 (however, the luminance decreases). Normally, γ = 0.1 is set.
(3) The 16-bit signal value v is converted into the 8-bit signal value v8 by the following calculation formula (compressed in 256 stages within the range from the lower limit value Lmin to the upper limit value Lmax).
v8 = (v−Lmin) · 255 / (Lmax−Lmin)
However, when v8> 255, it is saturated to v8 = 255, and when v8 <0, it is saturated to v8 = 0.
(4) The 16-bit color map Cmap (v, 3) (see (Equation 12)) is adjusted as follows.
When v8> Sv: Cmap (v, 3) ← Cmap (v, 3) (v8−Sv) ^{δ When} v8 ≦ Sv: No adjustment

  Cmap (v, n) (0 ≦ Cmap (v, n) ≦ 255, −32768 ≦ v ≦ 32767, n = 0 (R), 1 (G), 2 (B), 3 (α)) Yes, the 8-bit equivalent threshold is Sv, and the attenuation coefficient is δ (δ is a negative real number, usually δ = −0.5). The threshold value Sv is a threshold value of a signal value for attenuating the opacity, and is set to about Sv = 128 when the bone region is attenuated.

  On the other hand, when the signal value of the tomographic image group Do is 8 bits (see (Expression 2)), the control unit 11 adjusts the 8-bit color map Cmap (see (Expression 13)) as follows. .

(Formula 15)
When v> Sv: Cmap (v, 3) ← Cmap (v, 3) (v−Sv) ^δ
When v ≦ Sv: No adjustment

  Note that Cmap (v, n) (0 ≦ Cmap (v, n) ≦ 255, 0 ≦ v ≦ 255, n = 0 (R), 1 (G), 2 (B), 3 (α)). The threshold value is Sv, and the attenuation coefficient is δ (δ is a negative real number, usually δ = −0.5). The threshold value Sv is a threshold value of a signal value for attenuating the opacity, and is set to about Sv = 128 when the bone region is attenuated.

FIG. 22A is a diagram illustrating an attenuation characteristic of (v8−Sv) ^δ multiplied by Cmap (v, 3) (opacity curve) in (Expression 14). As shown in the figure, when the signal value v8> the threshold value Sv, the opacity is attenuated by the attenuation ratio (v8−Sv) ^δ corresponding to the difference value (v8−Sv) between the signal value v8 and the threshold value Sv. For example, as shown in FIG. 22B, when there is no difference in the distance from the center of the visceral region and the bone region, even if correction based on the correction table in step S51 described later is added to the opacity of both, The visceral region covers the bone region. Therefore, as shown in the figure, the opacity of the bone region is effectively reduced by using the difference between the signal values of the visceral region and the bone region and setting the threshold value Sv at the boundary between the signal values of the visceral region and the bone region. It can be attenuated and the visceral area is clearly visible.

  On the other hand, when the color map Cmap is not adjusted (step S49; “No”), the control unit 11 omits the above-described adjustment of the color map Cmap, and proceeds to step S51.

  The control unit 11 (voxel creation unit 31) obtains the tomographic image group Do acquired in step S41 (or the tomographic image group Do compressed in gradation in step S42) and the color map Cmap acquired in step S48 (or step S50). The voxel structure V (voxel data) is created based on the color map Cmap adjusted in step S43, the correction table S created in step S43, and the region of interest ROI specified in step S47 (step S51 in FIG. 21). ).

  In the present embodiment, the voxel structure V (voxel data) includes a voxel structure Vc that holds a 24-bit color value (RGB value) and a voxel structure Vα that holds an 8-bit opacity (α value). And two voxel structures. By creating the voxel structure in two as described above, the voxel processing relating to the color value and the opacity can be individually executed, so that the memory access and the like can be improved.

FIG. 23 is a flowchart showing a process for creating the voxel structure V (Vc, Vα).
First, based on the color map Cmap (v, n) (0 ≦ n ≦ 2) (= color palette) and the region of interest ROI, the control unit 11 stores a voxel structure Vc (x , Y, z, n) (0 ≦ V (x, y, z, n) ≦ 255, 0 ≦ n ≦ 2, 0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz -1; Resolution: Rxy, Rz) is created as follows (step S61 in FIG. 23).

(Formula 16)
(1) When Xs <x <Xe, Ys <y <Ye, and Zs <z <Ze are all satisfied (in the region of interest ROI)
Vc (x, y, z, n) = Cmap (Do (x, y, z), n) (0 ≦ n ≦ 2)
(2) When x ≦ Xs, x ≧ Xe, y ≦ Ys, y ≧ Ye, z ≦ Zs, or z ≧ Ze (when outside the region of interest ROI)
Vc (x, y, z, n) = 0 (0 ≦ n ≦ 2)

  Next, based on the color map Cmap (v, 3) (= opacity curve), the correction table S, and the region of interest ROI, the control unit 11 stores a voxel structure Vα (x, y) that holds 8-bit opacity. , Z) (0 ≦ Vα (x, y, z) ≦ 255, 0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1; resolution: Rxy, Rz) (Step S62 in FIG. 23).

(Formula 17)
(1) When Xs <x <Xe, Ys <y <Ye, and Zs <z <Ze are all satisfied (in the region of interest ROI)
Vα (x, y, z) = Cmap (Do (x, y, z), 3)
Vα (x, y, z) = Vα (x, y, z) · S (x, y, z)
(2) When x ≦ Xs, x ≧ Xe, y ≦ Ys, y ≧ Ye, z ≦ Zs, or z ≧ Ze (when outside the region of interest ROI)
Vα (x, y, z) = 0

  From the processing of Vα (x, y, z) = Vα (x, y, z) · S (x, y, z) in (Equation 17), the correction magnification S with respect to the opacity defined in the color map Cmap. Is multiplied to correct the opacity. In particular, the correction magnification S (x, y, z) decreases as the distance from the center of the ellipse approximating the outline of the subject region decreases, and thus the opacity of the voxel farther from the center of the subject region attenuates.

  For the purpose of improving memory access, the 32-bit voxel structure V is divided into a voxel structure Vc that holds color values and a voxel structure Vα that holds opacity. It does not have to be. In this case, the control unit 11 has one voxel structure V (x, y, z, n) (0 ≦ V (x, y, z, n) that holds a 32-bit color value and opacity as follows. n) ≦ 255, 0 ≦ n ≦ 3, 0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ z ≦ Sz−1; resolution: Rxy, Rz).

(Formula 18)
(1) When Xs <x <Xe, Ys <y <Ye, and Zs <z <Ze are all satisfied (in the region of interest ROI)
V (x, y, z, n) = Cmap (Do (x, y, z), n) (0 ≦ n ≦ 3)
V (x, y, z, 3) = V (x, y, z, 3) · S (x, y, z)
(2) When x ≦ Xs, x ≧ Xe, y ≦ Ys, y ≧ Ye, z ≦ Zs, or z ≧ Ze (when outside the region of interest ROI)
V (x, y, z, n) = 0 (0 ≦ n ≦ 3)

  The processing for creating the voxel structure shown in (Expression 16) to (Expression 18) described above is performed when the signal value of the tomographic image group Do is 16 bits (see (Expression 1)), and the color map Cmap corresponding to 16 bits is used. (Refer to (Expression 12)), and when the signal value of the tomographic image group Do is gradation-compressed to 8 bits (refer to (Expression 2)), a color map Cmap corresponding to 8 bits ((Expression 12)). 13)).

  After creating the voxel structure V (Vc, Vα), the control unit 11 performs a known smoothing process or a known shadow process on the voxel structure V (Vc, Vα) as necessary. May be.

  Subsequently, the control unit 11 (rendering unit 32) generates a volume rendering image based on the voxel structure V (step S52 in FIG. 21). Hereinafter, a rendering process by the ray casting method executed by the CPU (first rendering unit 32a) and a rendering process by the 3D texture mapping method executed by the GPU (second rendering unit 32b) will be described in order.

First, a rendering process by the ray casting method will be described with reference to FIGS.
FIG. 24 is a diagram showing an outline of the ray casting method. As shown in FIG. 24A, a ray is emitted from each point on the projection plane (rendered image) in an arbitrary direction with respect to the voxel data, and its trajectory from the deepest point (the deepest point where light reaches in the ray direction). Is reversed, and the luminance value of the point on the projection plane is calculated.
Also, as shown in FIG. 24 (b), in order to simplify the calculation, the voxel data is subjected to coordinate transformation, the ray is advanced from each point on the projection plane to the deepest point in the negative Z-axis direction, and the cumulative luminance calculation is simultaneously performed. Do.

  FIG. 25 is a flowchart showing the rendering process. First, the control unit 11 sets coordinate conversion parameters (step S71). A coordinate conversion process from the viewpoint coordinate system to the voxel coordinate system is performed in each process of a search control mask calculation (step S72) and a ray casting process (step S73) described later. Therefore, the control unit 11 sets and calculates in advance coordinate conversion parameters common to the coordinate conversion processing as follows.

<Coordinate transformation parameters>
Rotation parameter matrix R:
R = [R11 R12 R13;
R21 R22 R23;
R31 R32 R33]
(An inverse matrix of a 3 × 3 matrix for performing coordinate conversion from the voxel coordinate system to the viewpoint coordinate system. The GUI side instructs coordinate conversion from the voxel coordinate system to the viewpoint coordinate system. (Converts coordinates from to the voxel coordinate system.)
X-axis offset Xoff
Y-axis offset Yoff
Z-axis offset Zoff
ROI in the X-axis direction: Xs-Xe (0 ≦ Xs, Xe ≦ Sx−1)
ROI in the Y-axis direction: Ys−Ye (0 ≦ Ys, Ye ≦ Sy−1)
ROI in the Z-axis direction: Zs-Ze (0 ≦ Zs, Ze ≦ Sz−1)
Perspective transformation parameter: distance Dist in the Z-axis direction from the origin (positive real value, Dist = 0 for parallel projection)
Enlargement / reduction ratio Scale (same in XYZ axis direction)
Z-direction variable magnification Scz = Rxy / Rz
Z direction size after zooming in the Z direction Sz ′ = Sz · Scz
Coordinate transformation subsample offset: X-axis direction dx, Y-axis direction dy, Z-axis direction dz

The control unit 11 sets a unit matrix as an initial value for the rotation parameter matrix R described above. That is, R11 = 1, R12 = 0, R13 = 0, R21 = 0, R22 = 1, R23 = 0, R31 = 0, R32 = 0, and R33 = 1.
Then, according to the instruction of the GUI, any one of the X-axis center rotation Rx, the Y-axis center rotation Ry, and the Z-axis center rotation Rz (angle unit: radians) is sequentially specified, and a rotation matrix A is generated as follows. Then, the rotation parameter matrix R is multiplied from the right to update the rotation parameter matrix R. As a result, the inverse matrix of the rotation matrix from the voxel coordinate system to the viewpoint coordinate system generated by the GUI instruction is calculated.

Rotation matrix A
A = [A11 A12 A13;
A21 A22 A23;
A31 A32 A33]
Then,
Each element of the rotation matrix A in the case of the X-axis center rotation Rx is
A11 = 1, A12 = 0, A13 = 0
A21 = 0, A22 = cosRx, A23 = sinRx
A31 = 0, A32 = sinRx, A33 = cosRx
Each element of the rotation matrix A in the case of the Y-axis center rotation Ry is
A11 = cosRy, A12 = 0, A13 = sinRy
A21 = 0, A22 = 1, A23 = 0
A31 = −sinRy, A32 = 0, A33 = cosRy
Each element of the rotation matrix A in the case of the Z-axis center rotation Rz is
A11 = cosRz, A12 = sinRz, A13 = 0
A21 = −sinRz, A22 = cosRz, A23 = 0
A31 = 0, A32 = 0, A33 = 1
It becomes.
The rotation parameter matrix R is updated as R ← R × A.

  The coordinate conversion process described above is sequentially executed in each process of the search control mask calculation (step S72) and the ray casting process (step S73).

Next, the search control mask calculation process in step S72 of FIG. 25 will be described with reference to the flowchart of FIG.
The search control mask M (x, y, l) (0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1, 0 ≦ l ≦ L−1) is calculated for each pixel (x, y) of the rendering image to be calculated. Stores the Z coordinate of the leading opaque voxel. When performing sub-sampling (also referred to as anti-aliasing in order to remove jaggies in the volume rendering image) at the time of coordinate conversion, the control unit 11 performs the following processing L times (0 ≦ l ≦ L−1).

  First, the control unit 11 initializes the subsample offset value as dx = dy = dz = 0 by setting the number of times of processing to 1 = 0 (step S81).

  Subsequently, the control unit 11 performs the coordinate conversion from Zs = Sz′−1 to the pixel (x, y) having the values of x = 2m and y = 2n (m and n are integers). Search for an opaque voxel and record its Z coordinate in M (x, y, l). A similar search is performed for a pixel (x, y) having a value of either x = Sx-1 or y = Sy-1 (step S82). Thereby, as shown in FIG. 27A, the Z coordinate of the opaque voxel is calculated for the even-numbered pixels in the vertical and horizontal directions. The process for searching for opacity voxels will be described later with reference to FIG.

  Subsequently, the control unit 11 performs M (x−1, y) on a pixel (x, y) having values of x = 2m + 1 (x <Sx−1) and y = 2n (m and n are integers). , L) = M (x + 1, y, l), M (x, y, l) = M (x−1, y, l) is given, and M (x−1, y, l) ≠ M ( In the case of x + 1, y, l), the top opaque voxel is searched while performing coordinate conversion from Zs = Sz′−1 as described above (see FIG. 28), and the Z coordinate is recorded in M (x, y, l). (Step S83). Thereby, as shown in FIG. 27B, the Z coordinate of the opaque voxel is calculated for the odd-numbered pixels in the X direction. When the Z coordinates of the opaque voxels of the adjacent left and right pixels are the same, interpolation is performed using the Z coordinates (search for opaque voxels is not performed).

  Subsequently, the control unit 11 performs M (x, y−1) on all the pixels (x, y) having values of x (0 ≦ x ≦ Sx−1) and y = 2n + 1 (n is an integer). , L) = M (x, y + 1, l), M (x, y, l) = M (x, y−1, l) is given, and M (x, y−1, l) ≠ M ( In the case of x, y + 1, l), the top opaque voxel is searched while performing coordinate conversion from Zs = Sz′-1 as described above (see FIG. 28), and the Z coordinate is recorded in M (x, y, l). (Step S84). Thereby, as shown in FIG. 27C, the Z coordinate of the opaque voxel is calculated for the odd-numbered pixels in the Y direction. When the Z coordinate of the opaque voxel of the adjacent upper and lower pixels is the same, interpolation is performed with the Z coordinate (the search for the opaque voxel is not performed).

Subsequently, the control unit 11 updates the processing count to l ← l + 1, and updates the subsample offset values to dx ← dx + 1 / L, dy ← dy + 1 / L, and dz ← dz + 1 / L (step S85).
The control part 11 repeats the process of above-described step S82-S85 until it becomes the process frequency l = L-1 (step S86; Yes).

Next, the opacity voxel search process executed in the search control mask calculation process (FIG. 26) will be described with reference to the flowchart of FIG.
First, the control unit 11 initializes z to z = Zs and inputs a search target pixel (x, y) (step S91). Subsequently, the control unit 11 performs coordinate conversion on the three-dimensional coordinates (x, y, z), and calculates a voxel α value (n = 3) (step S92). Processing (coordinate conversion processing) for performing coordinate conversion and calculating the voxel α value will be described later.

  When the voxel α value calculated in step S92 is α = 0, the control unit 11 updates z ← z−m (for example, m = 8) (step S93). If z <0, the control unit 11 sets M (x, y, l) = − 1 (step S94) and ends the process. If z ≧ 0, the control unit 11 returns to step S92.

  On the other hand, when the voxel α value calculated in step S92 is α> 0, the control unit 11 first initializes i to i = 0 (step S95). Subsequently, the control unit 11 updates i ← i + 1 (step S96). When i <m and z + i <Zs, the control unit 11 performs coordinate conversion on the three-dimensional coordinates (x, y, z + i), and calculates the voxel α value. Is calculated (step S97). A process (coordinate conversion process) for performing the coordinate conversion and calculating the voxel α value will be described later.

  On the other hand, if i ≧ m or z + i ≧ Zs, the control unit 11 sets M (x, y, l) = z + i−1 (step S98), and ends the process.

  When the voxel α value calculated in step S97 is α = 0, the control unit 11 sets M (x, y, l) = z + i−1 (step S98) and ends the process.

  On the other hand, when the voxel α value calculated in step S97 is α> 0, the process returns to step S96, and the processes in steps S96 to S97 are repeated until the Z coordinate of the opaque voxel is set in step S98.

The coordinate conversion process executed in steps S92 and S97 of the opacity voxel search process (FIG. 28) will be described.
The coordinate conversion process is a process of converting the viewpoint coordinate system to the voxel coordinate system, and is the reverse of the conversion process on the GUI side. On the GUI side, it is assumed that clipping, scaling, Z-direction scaling processing, offset (simultaneously in the XYZ directions), rotation, and perspective transformation are performed in order of the ROI of interest, and the control unit 11 performs the three-dimensional view coordinate system given. The real number of the 8-bit voxel structure Vα (x, y, z) or 24-bit voxel structure Vc (x, y, z, n) corresponding to the coordinate value (x, y, z) (integer value) Coordinate values (xr, yr, zr) are calculated as follows.

The control unit 11 converts the coordinate value (x, y, z) (integer value) of the viewpoint coordinate system into a real value (xx, yy, zz) as follows.
xx = x-Sx / 2 + dx
yy = y-Sy / 2 + dy
zz = z−Sz ′ / 2 + dz

Subsequently, when the distance Dist> 0 in the Z-axis direction from the origin, the control unit 11 performs perspective transformation as follows.
zz ′ = Dist · zz / (Dist + zz)
xx ′ = xx · (Dist−zz ′) / Dist
yy '= yy. (Dist-zz') / Dist
Let (xx ′, yy ′, zz ′) after perspective transformation be (xx, yy, zz).

Subsequently, the control unit 11 performs a rotation process using the rotation parameter matrix R as follows.
xx ′ = R11 · xx + R12 · yy + R13 · zz
yy ′ = R21 · xx + R22 · yy + R23 · zz
zz ′ = R31 · xx + R32 · yy + R33 · zz
Let (xx ′, yy ′, zz ′) after the rotation process be (xx, yy, zz).

The control unit 11 simultaneously performs scaling, Z-direction scaling processing, and offset (simultaneously in the XYZ directions), and calculates the coordinate values (xr, yr, zr) (real values) of the corresponding voxels of the voxel structure as follows. get.
xr = xx / Scale + Sx / 2−Xoff
yr = yy / Scale + Sy / 2-Yoff
zr = zz / Scale / Scz + Sz / 2−Zoff

  Subsequently, the control unit 11 converts the coordinate value (xi, yi, zi) (integer value) by rounding off the decimal point to the calculated coordinate value (xr, yr, zr) (real value) of the voxel. And the fractions after the decimal point rounded down are (wx, wy, wz) (0 ≦ wx, wy, wz <1) (that is, xr = xi + wx, yr = yi + wy, zr = zi + wz).

  Then, the control unit 11 considers the ROI in the X-axis / Y-axis / Z-axis directions, and the voxel α value corresponding to the three-dimensional coordinate value (x, y, z) (integer value) of the viewpoint coordinate system ( Vα) is determined in three cases as follows.

1) When xi <Xs, xi> Xe, yi <Ys, yi> Ye, zi <Zs, or zi> Ze is satisfied (clipping range)
Vα = 0
2) When the above condition 1) is not satisfied, xi + 1> Xe, yi + 1> Ye, or zi + 1> Ze is satisfied (no interpolation)
Vα = Vα (xi, yi, zi)

3) When the above conditions 1) and 2) are not satisfied (interpolate)
Vα = (1-wz) (1-wy) (1-wx) · Vα (xi, yi, zi) + (1-wz) (1-wy) · wx · Vα (xi + 1, yi, zi) + ( 1-wz) * wy * (1-wx) * V [alpha] (xi, yi + 1, zi) + (1-wz) * wy * wx * V [alpha] (xi + 1, yi + 1, zi) + wz * (1-wy) (1- wx) · Vα (xi, yi, zi + 1) + wz · (1-wy) · wx · Vα (xi + 1, yi, zi + 1) + wz · wy · (1-wx) · Vα (xi, yi + 1, zi + 1) + wz · wy・ Wx ・ Vα (xi + 1, yi + 1, zi + 1)

  Note that when the RGB value of the voxel is acquired as in a ray casting process (see FIG. 29) described later, the control unit 11 sets the three-dimensional coordinate value (x, y, z) (integer value) of the viewpoint coordinate system. The RGB values (Vc (n)) (0 ≦ n ≦ 2) of the corresponding voxels are determined in three cases as follows.

1) When xi <Xs, xi> Xe, yi <Ys, yi> Ye, zi <Zs, or zi> Ze is satisfied (clipping range)
Vc (n) = 0 (0 ≦ n ≦ 2)
2) When the above condition 1) is not satisfied, xi + 1> Xe, yi + 1> Ye, or zi + 1> Ze is satisfied (no interpolation)
Vc (n) = Vc (xi, yi, zi, n) (0 ≦ n ≦ 2)

3) When the above conditions 1) and 2) are not satisfied (interpolate)
Vc (n) = (1-wz) (1-wy) (1-wx) .Vc (xi, yi, zi, n) + (1-wz) (1-wy) .wx.Vc (xi + 1, yi , Zi, n) + (1-wz) · wy · (1-wx) · Vc (xi, yi + 1, zi, n) + (1-wz) · wy · wx · Vc (xi + 1, yi + 1, zi, n ) + Wz · (1−wy) (1−wx) · Vc (xi, yi, zi + 1, n) + wz · (1−wy) · wx · Vc (xi + 1, yi, zi + 1, n) + wz · wy · (1 −wx) · Vc (xi, yi + 1, zi + 1, n) + wz · wy · wx · Vc (xi + 1, yi + 1, zi + 1, n) (0 ≦ n ≦ 2)

Next, the ray casting process in step S73 of FIG. 25 will be described with reference to the flowchart of FIG.
First, the control unit 11 sets all initial values of the generated 24-bit (RGB) rendered image Image (x, y, n) to 0 (Image (x, y, n) = 0, n = 0 ( R), 1 (G), 2 (B)). Then, the following processing is executed for each two-dimensional coordinate (x, y) (0 ≦ x ≦ Sx−1, 0 ≦ y ≦ Sy−1) as the number of subsamples L.

  First, the control unit 11 sets the light source Light (n) (n = 0 (R), 1 (G), 2 (B), 0 ≦ Light (n) ≦ 255), and sets the subsample offset value to dx. = Dy = dz = 0 and the number of processing times 1 = 0 is initialized (step S101).

  Subsequently, the control unit 11 sets z = M (x, y, l) (step S102). When z <0, the control unit 11 proceeds to step S110 and determines the RGB pixel value at the pixel (x, y).

  On the other hand, if z ≧ 0, the control unit 11 initializes the virtual light intensity Trans = 1.0 and the accumulated luminance value Energy (n) = 0.0 (0 ≦ n ≦ 2) (step S103) and three-dimensional The coordinates (x, y, z) are transformed to obtain the voxel α value (Vα), and Alpha = Vα / 255 is calculated (step S104).

  When Alpha calculated in step S104 is Alpha = 0 and z = 0, the process proceeds to step S110, and the RGB pixel value in the pixel (x, y) is determined.

  On the other hand, when the Alpha calculated in Step S104 is Alpha> 0, the control unit 11 performs coordinate conversion on the three-dimensional coordinates (x, y, z) to obtain the voxel RGB value Vc (n) (Step S107). The RGB value Vc (n) of the voxel is acquired by the coordinate conversion process described above.

  Subsequently, the control unit 11 updates the accumulated luminance as Energy (n) = Energy (n) + Trans / Alpha · Vc (n) / 255 and the transmitted light intensity as Trans = Trans · (1.0−Alpha) ( Step S108). When Alpha = 1.0 or Trans <0.001, the process proceeds to step S110, and the RGB pixel value at the pixel (x, y) is determined. Otherwise, the control unit 11 updates z = z−1 (step S109). If z ≧ 0, the control unit 11 returns to step S104. If z ≦ 0, the control unit 11 proceeds to step S110, and the pixel (x, Determine RGB pixel values in y).

  On the other hand, when the Alpha calculated in Step S104 is Alpha = 0 and z> 0, the control unit 11 performs coordinate conversion from Zs = z−1 and searches for the Z coordinate z ′ of the leading opaque voxel (Step S105). . The search for opaque voxels is executed by the above-described opaque voxel search process (see FIG. 28).

  When the searched z ′ is z ′ <0, the process proceeds to step S110, and the RGB pixel value at the pixel (x, y) is determined. When z ′ ≧ 0, the control unit 11 updates z = z ′ (step S106) and returns to the process of step S104.

In step S110, the control unit 11 determines the RGB pixel value in the pixel (x, y) as follows.
Image (x, y, n) = Image (x, y, n) + k · Energy (n) · Light (n) / L
Here, k is an intensity magnification, and an initial value is set to k = 1.0.

  Subsequently, the control unit 11 updates the number of processes to l ← l + 1, and updates the subsample offset values to dx ← dx + 1 / L, dy ← dy + 1 / L, and dz ← dz + 1 / L (step S111). The control unit 11 repeats the processing of steps S102 to S111 described above until the number of processing times l = L−1 (step S112; Yes), and the final rendered image Image (x, y, n) (n ≦ 0). ≦ 2) is obtained.

  The rendering process by the ray casting method executed by the CPU (first rendering unit 32a) has been described above with reference to FIGS.

  Next, a rendering process by the 3D texture mapping method executed by the GPU (second rendering unit 32b) will be described with reference to FIGS.

  Since rendering by the GPU is performed using the OpenGL system, the control unit 11 registers the created voxel structure V as a 3D texture in the OpenGL system. Specifically, the control unit 11 converts the 24-bit voxel structure Vc (see (Expression 16)) and the 8-bit voxel structure Vα (see Expression (17)) into a 32-bit voxel structure V ((Expression 16). 18), 3D texture is generated using the glTexImage3D function, and the 3D texture is registered in the OpenGL system.

  In addition, as shown in FIG. 30, the control unit 11 defines each slice constituting the 3D texture as a rectangle, sets a stacked rectangle arranged in the Z-axis direction on the XY coordinate plane of the three-dimensional space, and sets the texture coordinates. And 3D texture are pasted on each square. In the 3D texture mapping method, a rendering image in an arbitrary direction is generated by performing coordinate conversion on a 3D texture map attached to a quadrangle in accordance with the line-of-sight direction.

  In the rendering process by the GPU (second rendering unit 32b), moire may occur during rendering. In order to suppress the occurrence of moire during rendering, the control unit 11 performs in advance color correction processing of transparent voxels on the ROI boundary surface of the region of interest as follows.

  FIG. 31 is a flowchart of the color correction process. The control unit 11 extracts boundary surface voxels (x, y, z) (step S121). The boundary surface voxel (x, y, z) is a voxel satisfying any of x = Xs, x = Xe, y = Ys, y = Ye, z = Zs, or z = Ze.

  Subsequently, the control unit 11 extracts opaque voxels from 27 neighboring voxels of the extracted boundary surface voxels, and calculates the average colors (Ra, Ga, Ba) of the extracted opaque voxels (step S122). An opaque voxel is a voxel that satisfies the condition of Vα (x + i, y + j, z + k)> 0 (−1 ≦ i ≦ 1, −1 ≦ j ≦ 1, −1 ≦ k ≦ 1). Note that, when an opaque voxel is not extracted, the control unit 11 sets Ra = Ga = Ba = 0.

Then, the control unit 11 sets the opacity of the boundary surface voxel to 0 and sets the RGB value to the average color of the opaque voxels as described below (step S123).
Vα (x, y, z) = 0
Vc (x, y, z, 1) = Ra, Vc (x, y, z, 2) = Ga, Vc (x, y, z, 3) = Ba

  The control unit 11 executes the color correction process (steps S121 to S123) described above for all the boundary surface voxels (x, y, z).

  Next, the flow of rendering processing by the 3D texture mapping method will be described with reference to the flowchart of FIG. First, the control unit 11 sets an OpenGL perspective transformation parameter (step S131). When rendering is performed in the virtual endoscope mode, the control unit 11 sets the following perspective transformation parameters for the OpenGL system.

<Transparent transformation parameters>
-Camera viewing angle (unit: degrees, 0 for parallel projection, equivalent to focal length)
-Viewpoint coordinates (camera position, usually set on the Z axis)
・ Gaze point coordinates (fixed point of the subject, usually fixed at z = 0 on the Z axis)
・ Near clipping position (distance from viewpoint, slightly distance from viewpoint)
-Far clipping position (distance from viewpoint, usually fixed at infinity such as 1000 and not clipped)

  On the other hand, when normal rendering is performed, the viewing angle of the camera is set to 0 and the parallel projection mode is set.

Subsequently, the control unit 11 performs OpenGL projection screen setting (step S132).
FIG. 33 is a flowchart showing a projection screen setting process. As shown in FIG. 33, the control unit 11 sets the screen size (number of vertical and horizontal pixels, vertical and horizontal aspect ratio) of the rendered image (step S141), and performs parallel projection (normal appearance rendering) or perspective projection (endoscopic mode). ) Is set (step S142). When the perspective projection is set, the control unit 11 sets the perspective projection parameters (camera viewing angle (focal length), viewpoint coordinates, gazing point coordinates, near clipping position, far clipping position) ( Step S143).

Returning to the flowchart of FIG.
Subsequently, the control unit 11 sets OpenGL coordinate conversion parameters (rotation, scale, movement, Z-direction scaling parameter Rxy / Rz) as follows (step S133).

<OpenGL coordinate conversion parameters>
Rotation parameter matrix R:
R = [R11 R12 R13 R14;
R21 R22 R23 R24;
R31 R32 R33 R34;
R41 R42 R43 R44]
(This is a 4 × 4 matrix for coordinate conversion from the voxel coordinate system to the viewpoint coordinate system, and the rendering side also instructs OpenGL to perform coordinate conversion in the same direction as the GUI side. Actually, it is defined by a 3 × 3 coordinate conversion matrix. Yes, but extended to 4 × 4 coordinate transformation matrix due to OpenGL specifications)
X-axis offset Xoff
Y-axis offset Yoff
Z-axis offset Zoff
Enlargement / reduction ratio Scale (same in XYZ axis direction)
Z-direction variable magnification Scz = Rxy / Rz
Z direction size after zooming in the Z direction Sz ′ = Sz · Scz

The control unit 11 sets a unit matrix as an initial value for the rotation parameter matrix R described above. That is, R11 = 1, R12 = 0, R13 = 0, R14 = 0, R21 = 0, R22 = 1, R23 = 0, R24 = 0, R31 = 0, R32 = 0, R33 = 1, R41 = 0, R42 = 0, R43 = 0, and R44 = 1 are set.
Then, according to the instruction of the GUI, any one of the X-axis center rotation Rx, the Y-axis center rotation Ry, and the Z-axis center rotation Rz (angle unit: radians) is sequentially specified, and a rotation matrix A is generated as follows. Then, the rotation parameter matrix R is multiplied from the left to update the rotation parameter matrix R.

Rotation matrix A
A = [A11 A12 A13 A14;
A21 A22 A23 A24;
A31 A32 A33 A34;
A41 A42 A43 A44]
Then,
Each element of the rotation matrix A in the case of the X-axis center rotation Rx is
A11 = 1, A12 = 0, A13 = 0, A14 = 0
A21 = 0, A22 = cosRx, A23 = sinRx, A24 = 0
A31 = 0, A32 = sinRx, A33 = cosRx, A34 = 0
A41 = 0, A42 = 0, A43 = 0, A44 = 1
Each element of the rotation matrix A in the case of the Y-axis center rotation Ry is
A11 = cosRy, A12 = 0, A13 = sinRy, A14 = 0
A21 = 0, A22 = 1, A23 = 0, A24 = 0
A31 = −sinRy, A32 = 0, A33 = cosRy, A34 = 0
A41 = 0, A42 = 0, A43 = 0, A44 = 1
Each element of the rotation matrix A in the case of the Z-axis center rotation Rz is
A11 = cosRz, A12 = sinRz, A13 = 0, A14 = 0
A21 = −sinRz, A22 = cosRz, A23 = 0, A24 = 0
A31 = 0, A32 = 0, A33 = 1, A34 = 0
A41 = 0, A42 = 0, A43 = 0, A44 = 1
It becomes.
The rotation parameter matrix R is updated as R ← A × R.

Returning to the flowchart of FIG.
The control unit 11 executes rendering processing using OpenGL (step S134).
FIG. 34 is a flowchart showing rendering processing by OpenGL. As shown in FIG. 34, the control unit 11 performs scaling and Z-direction scaling processing on the 3D texture map based on the enlargement / reduction ratio Scale (same in the XYZ-axis directions) and the Z-direction scaling ratio Scz ( In step S151), a rotation process is performed on the 3D texture map based on the rotation parameter matrix R (step S152), and based on the offset Xoff in the X-axis direction, the offset Yoff in the Y-axis direction, and the offset Zoff in the Z-axis direction. By performing offset processing on the 3D texture map (step S153), predetermined coordinate transformation is performed on the 3D texture map.

  Subsequently, the control unit 11 sets a stacked quadrangle in which quadrangles on the XY coordinate plane of the three-dimensional space are arranged in the Z-axis direction, and associates the 3D texture map (see step S154 and FIG. 30). Then, the control unit 11 performs rendering processing using OpenGL, performs scan conversion while pasting the 3D texture map on the stacked quadrangle, and performs alpha blending processing from the back to the Z direction (viewpoint direction) on the frame memory (step S155). ). That is, the color value (RGB value) of the voxel of the 3D texture map after coordinate conversion corresponding to the rectangular XY coordinates on the line of sight parallel to the Z-axis direction from a predetermined viewpoint is based on the opacity (α value) of the voxel. Are obtained by alpha blending in the order of a rectangle far from the viewpoint and given as a pixel value of a volume rendering image.

  The rendering process by the 3D texture mapping method executed by the GPU (second rendering unit 32b) has been described above with reference to FIGS.

  The control unit 11 performs a rendering process based on the ray casting method (FIG. 25) executed by the CPU (first rendering unit 32a) or a rendering process based on the 3D texture mapping method (FIG. 25) executed by the GPU (second rendering unit 32b). 32), a volume rendering image is generated, and the generated volume rendering image is displayed on the display unit 16.

Returning to the flowchart of FIG.
The control unit 11 (data output unit 26) can output various data and store it in the storage unit 12 by receiving a data storage operation from the user (step S53 in FIG. 21). For example, the control unit 11 outputs and stores the MPR image generated and displayed in step S44 and the volume rendering image (three-dimensional image) generated and displayed in step S52.

[Example]
Finally, the rendered images generated by the conventional method and the proposed method are compared. The conventional method is a method of creating a voxel structure V (Vc, Vα) by applying the color map Cmap acquired in step S48 of FIG. 21 as it is to the tomographic image group Do and generating a rendering image. That is, in the conventional method, the voxel structure V (Vc, Vα) is created based on the following (Equation 19).

(Formula 19)
(1) When Xs <x <Xe, Ys <y <Ye, and Zs <z <Ze are all satisfied (in the region of interest ROI)
Vc (x, y, z, n) = Cmap (Do (x, y, z), n) (0 ≦ n ≦ 2)
Vα (x, y, z) = Cmap (Do (x, y, z), 3)
(2) When x ≦ Xs, x ≧ Xe, y ≦ Ys, y ≧ Ye, z ≦ Zs, or z ≧ Ze (when outside the region of interest ROI)
Vc (x, y, z, n) = 0 (0 ≦ n ≦ 2)
Vα (x, y, z) = 0

The proposed method is a method of generating a rendering image based on the present embodiment (a method of generating a rendering image by the processing of FIG. 21). Here, parameters for calculating the correction magnification S (see (Equation 9) (Equation 10)) are k = 1.1 (attenuation coefficient in the XY plane), kz = 0.0 (attenuation coefficient in the Z direction). , Amp = ampz = 1.0 (amplification magnification), xoffset = 2 (offset in the X direction), yoffset = 5 (offset in the Y direction), xoffset = 0 (offset in the Z direction), a = 0.58 (horizontal) Direction size correction coefficient), b = 0.58 (vertical size correction coefficient), and c = 1.0 (slice total magnification).
In addition, neither the proposed method nor the conventional method performs the adjustment of the color map Cmap (the process of step S50 in FIG. 21), and both perform rendering using the first rendering unit 32a (second rendering unit 32b). Can be used to obtain almost the same rendered image).

  FIG. 35 is an image obtained by observing the rendered image generated by the conventional method from the front, and FIG. 36 is an image obtained by observing the rendered image generated by the proposed method from the front. In FIG. 35 (conventional technique), the ribs and sternum are not sufficiently transparent, and it is difficult to observe internal organs and the like, whereas in FIG. The whole is transparent, and internal organs and the like can be observed well. FIG. 37 and FIG. 38 are images obtained by observing the rendered images generated by the proposed method from the back and side surfaces, respectively. As shown in FIGS. 37 and 38, although the vertebra remains, it can be seen that the rib is satisfactorily transparent as in FIG.

  The fourth embodiment has been described above. According to the fourth embodiment, a volume rendering image is generated and displayed based on the tomographic image group Do. In particular, when the voxel structure V (voxel data) is created by converting the signal value v of each pixel (x, y, z) of the tomographic image group Do into color values and opacity with reference to the color map Cmap. Rather than using the opacity defined in the color map Cmap (v, 3) (opacity curve) as it is, correction calculated based on ellipse information (ellipse parameter P (z)) approximating the outline of the subject area The opacity is corrected by multiplying by the magnification S (x, y, z) (see (Equation 17)). Then, a voxel structure (voxel data) is created using the corrected opacity, and a rendering image is generated. Thereby, it is possible to obtain a rendering image in which a specific area is clearly visualized without creating a three-dimensional mask as in the conventional case. In this disclosure, it is necessary to create a correction table for the correction magnification S (x, y, z) in place of the conventional three-dimensional mask. However, these methods do not take a method of creating a correction table in advance as described above. In addition, a method of calculating each of the voxels by (Equation 9) and (Equation 10) can be taken at the time of voxel creation, and it is not necessary to secure them in the memory like a three-dimensional mask.

  Further, the bone region can be attenuated more effectively by adjusting the color map Cmap (v, 3) (opacity curve) according to the signal value (see FIG. 21).

  In addition, as in conventional mask processing, in the visible / invisible control (on / off control) of each voxel, a boundary region in which a visible internal region and an invisible bone region are mixed in a single voxel is expressed. In some cases, a natural rendered image could not be obtained. On the other hand, in the present embodiment, the opacity is gradually attenuated toward the outside of the subject area instead of the on / off control as in the mask process, and a smooth state is obtained while providing an intermediate state between on and off. By performing the control, a natural rendered image can be obtained.

  The preferred embodiments of the medical image processing apparatus and the like according to the present disclosure have been described above with reference to the accompanying drawings, but the present disclosure is not limited to such examples. It will be apparent to those skilled in the art that various changes and modifications can be made within the scope of the technical idea disclosed in the present application, and these naturally belong to the technical scope of the present disclosure. Understood.

1, 1a, 1b: medical image processing device 21: image acquisition unit 22: gradation compression unit 23: correction table creation unit 23a: geometric information acquisition unit 23b: geometric information correction unit 23c: correction magnification calculation unit 24: MPR generation unit 25: MPR display unit 26: data output unit 27: parameter adjustment unit 28: area designation unit 29: color map acquisition unit 30: color map adjustment unit 31: voxel creation unit 32: rendering unit 32a: first rendering unit 32b: first 2 Rendering unit 40: Parameter adjustment screen S: Correction magnification, correction table

Claims (26)

Image acquisition means for acquiring a plurality of tomographic images;
Geometric information acquisition means for approximating the outline of the subject area of the tomographic image with a predetermined geometric shape for each tomographic image, and acquiring information on the geometric shape;
MPR image generation means for correcting a signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and generating an MPR image based on the tomographic image whose signal value is corrected; A medical image processing apparatus comprising: Correction magnification calculating means for calculating, for each tomographic image, a correction magnification corresponding to each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image;
The medical image processing apparatus according to claim 1, wherein the MPR image generation unit multiplies the signal value of each pixel of the tomographic image by the correction magnification calculated by the correction magnification calculation unit. The correction magnification calculation means further comprises table creation means for creating a correction table for storing the calculated correction magnification,
The medical image processing apparatus according to claim 2, wherein the MPR generation unit multiplies a signal value of each pixel of the tomographic image by a correction magnification obtained by referring to the correction table. The medical image processing apparatus according to claim 2, wherein the correction magnification calculation unit decreases the correction magnification as the distance from the geometric center of the geometric shape increases. The geometric shape is an ellipse;
The geometric information acquisition means acquires the center coordinates of the ellipse, the horizontal size of the ellipse, and the vertical size of the ellipse as the geometric shape information,
5. The medical image processing apparatus according to claim 2, wherein the correction magnification calculation unit calculates the correction magnification based on the center coordinate, the horizontal size, and the vertical size. 6. . The correction magnification calculation means maximizes the correction magnification in the center coordinate, and the correction on each concentric ellipse inside the ellipse in which the ratio of the central coordinate to the horizontal size and the vertical size is the same. The medical image processing apparatus according to claim 5, wherein the magnification is reduced as the horizontal size and the vertical size of each concentric ellipse increase. The correction magnification calculation means corrects the center coordinates by adding a predetermined offset to the center coordinates, and the horizontal size and the vertical size to the predetermined horizontal size and the vertical size. The horizontal size and the vertical size are corrected by multiplying by the magnification with respect to, and the correction magnification is calculated based on the corrected center coordinates, the horizontal size, and the vertical size. The medical image processing apparatus according to claim 5 or 6. The medical image processing apparatus according to claim 7, further comprising an adjustment unit that allows a user to adjust the offset, the magnification with respect to the size in the horizontal direction, and the magnification with respect to the size in the vertical direction. 2. The medical image according to claim 1, wherein the MPR generation unit corrects the signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and information on a slice position of the tomographic image. Processing equipment. Correction magnification calculating means for calculating, for each tomographic image, a correction magnification corresponding to each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and information on a slice position of the tomographic image; In addition,
The medical image processing apparatus according to claim 9, wherein the MPR image generation unit multiplies the signal value of each pixel of the tomographic image by the correction magnification calculated by the correction magnification calculation unit. The correction magnification calculation means further comprises table creation means for creating a correction table for storing the calculated correction magnification,
The medical image processing apparatus according to claim 10, wherein the MPR generation unit multiplies a signal value of each pixel of the tomographic image by a correction magnification obtained by referring to the correction table. 12. The medical image according to claim 10, wherein the correction magnification calculation unit decreases the correction magnification as the distance from the geometric center of the geometric shape decreases, and decreases as the slice position is located closer to the end from the center. Processing equipment. The geometric shape is an ellipse;
The geometric information acquisition means acquires the center coordinates of the ellipse, the horizontal size of the ellipse, and the vertical size of the ellipse as the geometric shape information,
The information of the slice position is the total number of slices and the slice order,
The correction magnification calculation unit calculates the correction magnification based on the center coordinates, the horizontal size, the vertical size, the total number of slices, and the slice order. The medical image processing apparatus described in 1. The correction magnification calculation means maximizes the correction magnification in the center coordinate, and the correction on each concentric ellipse inside the ellipse in which the ratio of the central coordinate to the horizontal size and the vertical size is the same. The magnification is decreased as the horizontal size and vertical size of each concentric ellipse are increased, and the correction magnification is decreased as the difference between the total number of slices / 2 and the slice rank is increased. The medical image processing apparatus described. The correction magnification calculation means corrects the center coordinates by adding a predetermined offset to the center coordinates, and the horizontal size and the vertical size to the predetermined horizontal size and the vertical size. The horizontal size and the vertical size are corrected by multiplying by the magnification with respect to, and the total number of slices is corrected by multiplying by the magnification with respect to the predetermined total number of slices, and the slice order is added with a predetermined slice offset The correction magnification is calculated based on the corrected center coordinates, the horizontal size, the vertical size, the total number of slices, and the slice order. Medical image processing apparatus. The medical image processing according to claim 15, further comprising: an adjustment unit that allows a user to adjust the offset, a magnification with respect to the horizontal size, a magnification with respect to the vertical size, a magnification with respect to the total number of slices, and the slice offset. apparatus. The approximate accuracy of the geometric shape of the tomographic image is determined based on a predetermined criterion, and the geometric shape information of the tomographic image that does not satisfy the predetermined criterion is based on the geometric shape information of the other tomographic image that satisfies the predetermined criterion. The medical image processing apparatus according to claim 1, further comprising: a geometric information correcting unit that corrects the information. The geometric shape is an ellipse;
The geometric information acquisition means acquires the center coordinates of the ellipse, the horizontal size of the ellipse, and the vertical size of the ellipse as the geometric shape information,
The geometric information correcting unit is configured to reduce the horizontal size of the ellipse in the tomographic image in which the horizontal size of the ellipse matches the horizontal size of the image, and the horizontal size of the ellipse is smaller than the horizontal size of the image. The medical image processing apparatus according to claim 17, wherein correction is performed so as to be equal to a ratio of a horizontal size and a vertical size of an ellipse in the other tomographic image. Color map acquisition means for acquiring a color map that defines the correspondence between signal values, color values, and opacity;
By referring to the color map, the signal value of each pixel of the tomographic image is converted into a color value and opacity corresponding to the signal value, and voxel data that is a set of voxels that retain the color value and opacity is obtained. Voxel creation means to create,
The medical image processing apparatus according to claim 1, further comprising: a rendering unit that generates a volume rendering image based on the voxel data. The medical image processing apparatus according to claim 19, wherein the voxel creation unit corrects the opacity according to a signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image. 21. The medical image processing apparatus according to claim 20, wherein the voxel creation unit further corrects the opacity according to a signal value of each pixel of the tomographic image based on information on a slice position of the tomographic image. Color map adjusting means for adjusting the color map so as to attenuate the opacity corresponding to the signal value according to a difference value between the signal value and a predetermined threshold;
The medical image processing apparatus according to any one of claims 19 to 21, wherein the voxel creation unit creates voxel data with reference to the adjusted color map. Color map creating means for creating a color map that defines the correspondence between the signal value of a predetermined range, the color value, and the opacity based on the correspondence between the representative signal value, the color value, and the opacity. Prepared,
The medical image processing apparatus according to any one of claims 19 to 22, wherein the color map acquisition unit acquires the created color map. The computer controller
An image acquisition step of acquiring a plurality of tomographic images;
For each tomographic image, a geometric information acquisition step of approximating the outline of the subject area of the tomographic image with a predetermined geometric shape, and acquiring information on the geometric shape;
An MPR image generation step of correcting a signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and generating an MPR image based on the tomographic image in which the signal value is corrected; A medical image processing method. Computer
Image acquisition means for acquiring a plurality of tomographic images;
Geometric information acquisition means for approximating the outline of the subject area of the tomographic image with a predetermined geometric shape for each tomographic image, and acquiring information on the geometric shape;
An MPR image generating unit that corrects the signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image, and generates an MPR image based on the tomographic image in which the signal value is corrected;
Program to function as. An MPR image generation method executed by an MPR image generation device that generates an MPR image based on a tomographic image,
An image acquisition step of acquiring a plurality of tomographic images;
For each tomographic image, a geometric information acquisition step of approximating the outline of the subject area of the tomographic image with a predetermined geometric shape, and acquiring information on the geometric shape;
An MPR image generation step of correcting a signal value of each pixel of the tomographic image based on information on a geometric shape corresponding to the tomographic image and generating an MPR image based on the tomographic image in which the signal value is corrected;
An MPR image generation method including: