CN116977391A - Underwater three-dimensional reconstruction method based on binocular multi-line structured light - Google Patents
Underwater three-dimensional reconstruction method based on binocular multi-line structured light Download PDFInfo
- Publication number
- CN116977391A CN116977391A CN202310948944.0A CN202310948944A CN116977391A CN 116977391 A CN116977391 A CN 116977391A CN 202310948944 A CN202310948944 A CN 202310948944A CN 116977391 A CN116977391 A CN 116977391A
- Authority
- CN
- China
- Prior art keywords
- camera
- underwater
- point cloud
- imu
- image
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 45
- 238000012545 processing Methods 0.000 claims abstract description 26
- 230000001360 synchronised effect Effects 0.000 claims abstract description 10
- 230000010354 integration Effects 0.000 claims abstract description 8
- 239000011159 matrix material Substances 0.000 claims description 49
- 238000005259 measurement Methods 0.000 claims description 24
- 238000003384 imaging method Methods 0.000 claims description 21
- 239000013598 vector Substances 0.000 claims description 21
- 230000005484 gravity Effects 0.000 claims description 19
- 238000013519 translation Methods 0.000 claims description 16
- 238000012937 correction Methods 0.000 claims description 14
- 238000001914 filtration Methods 0.000 claims description 7
- 230000009466 transformation Effects 0.000 claims description 5
- 238000009499 grossing Methods 0.000 claims description 4
- 238000007781 pre-processing Methods 0.000 claims description 4
- 230000009467 reduction Effects 0.000 claims description 4
- 238000012216 screening Methods 0.000 claims description 4
- 230000008569 process Effects 0.000 abstract description 7
- 230000003287 optical effect Effects 0.000 abstract description 4
- 230000007547 defect Effects 0.000 abstract description 3
- 238000012360 testing method Methods 0.000 abstract description 2
- 238000001514 detection method Methods 0.000 description 7
- 238000005516 engineering process Methods 0.000 description 7
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 4
- 230000001133 acceleration Effects 0.000 description 3
- 230000008859 change Effects 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000002708 enhancing effect Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000001131 transforming effect Effects 0.000 description 2
- 230000003044 adaptive effect Effects 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 239000000284 extract Substances 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000012634 optical imaging Methods 0.000 description 1
- 239000002245 particle Substances 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000013535 sea water Substances 0.000 description 1
- 238000012876 topography Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/50—Depth or shape recovery
- G06T7/521—Depth or shape recovery from laser ranging, e.g. using interferometry; from the projection of structured light
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/10—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
- G01C21/12—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
- G01C21/16—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/50—Depth or shape recovery
- G06T7/55—Depth or shape recovery from multiple images
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10028—Range image; Depth image; 3D point clouds
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Theoretical Computer Science (AREA)
- Optics & Photonics (AREA)
- Automation & Control Theory (AREA)
- Length Measuring Devices By Optical Means (AREA)
- Image Processing (AREA)
Abstract
The invention belongs to the technical field of underwater environment three-dimensional sweep testing of underwater robots, and particularly relates to an underwater three-dimensional reconstruction method based on binocular multi-line structured light. The invention uses multi-line structured light to scan the underwater target with high precision, and effectively solves the defects of point cloud distortion, information loss and the like of three-dimensional reconstruction of lattice laser and single-line laser. The underwater target is scanned by the binocular camera with the multi-line laser transmitters and the narrow-band optical filters, two groups of photos are not required to be shot underwater for processing, and scanning efficiency is improved. According to the invention, through the binocular multi-line underwater three-dimensional reconstruction system arranged on the underwater robot, the IMU is utilized to perform pre-integration, the object to be measured is scanned in the moving process, each frame of photo of the underwater binocular camera comprises a plurality of lasers, so that the point cloud density is improved, the scanning frequency is reduced, and the synchronous three-dimensional reconstruction of the underwater environment is realized.
Description
Technical Field
The invention belongs to the technical field of underwater environment three-dimensional sweep testing of underwater robots, and particularly relates to an underwater three-dimensional reconstruction method based on binocular multi-line structured light.
Background
With the continuous development of ocean resource development, submarine topography detection and underwater scientific research, the underwater vision three-dimensional reconstruction technology plays an important role in the underwater target detection technology. The underwater vision three-dimensional reconstruction integrates an optical imaging technology, an image processing technology, a communication signal processing technology, a synchronous image construction and positioning technology and the like, and provides powerful guarantee for underwater target detection.
The current underwater detection technology mainly comprises sonar detection, wherein when sound waves propagate in water, the sound waves are influenced by the temperature, water pressure, salinity, suspended particles in the water and other environmental factors of the sea water; compared with sound waves, the laser has the advantages of high brightness, shorter pulse, higher collimation degree and the like, so that the detection distance measurement is higher, and the positioning and imaging are more accurate.
The traditional lattice or single-line laser scanning is easily affected by scattering and diffusion in the underwater propagation process, so that the problems of information loss, sparse point cloud data, low scanning speed, small detection range and the like are caused.
Disclosure of Invention
The invention aims to solve the problems of great influence of water quality on the current acoustics, information loss of underwater scanning of lattice laser and single-line structured light, sparse point cloud and the like, and provides an underwater three-dimensional reconstruction method based on binocular multi-line structured light.
An underwater three-dimensional reconstruction method based on binocular multi-line structured light comprises the following steps:
step 1: a binocular camera and a multi-line laser transmitter are arranged at the bottom of the underwater robot; the binocular camera comprises a left camera and a right camera which are respectively arranged at two sides of the multi-line laser transmitter;
step 2: calibrating the binocular camera and an inertial measurement unit IMU;
step 3: calibrating light planes of the left camera and the right camera;
step 4: aligning the binocular camera with an inertial measurement unit IMU time stamp, and carrying out synchronous data acquisition at the time t; pre-integrating an Inertial Measurement Unit (IMU), processing a laser image of a binocular camera, and extracting point cloud information;
step 5: the method comprises the steps that a left camera and a right camera are used for respectively obtaining a left image and a right image of an object to be measured, and the left image and the right image are subjected to three-dimensional correction, so that corrected left image and corrected right image are aligned; matching the corrected left image and right image to obtain a line laser imaging matching point pair; obtaining left and right image parallax according to the line laser imaging matching point pairs, and calculating three-dimensional data of each laser line according to the left and right image parallax;
step 6: performing pose matching correction based on multi-frame registration of the three-dimensional point cloud data, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame registration;
step 7: and processing the redundant three-dimensional data, and removing noise points and outliers by adopting a voxel filtering method to reconstruct the space three-dimensionally.
Further, calibrating the binocular camera in the step 2 includes calibrating a lens focal length f, a base line distance B, and an internal reference K of the left camera l Internal reference K of right camera r Distortion parameter [ k ] 1 k 2 p 1 p 2 k 3 ]A rotation matrix R and a translation matrix T of the binocular camera; wherein k is 1 、k 2 、k 3 For radial distortion coefficient, p 1 、p 2 Is a tangential distortion coefficient;
the calibrating the inertial measurement unit IMU comprises setting an IMU coordinate system as a world coordinate system, and then converting the binocular camera coordinate system into the IMU coordinate system into the following relationship:
wherein ,and->Is two-dimensional coordinates under binocular camera coordinates, +.>Is a three-dimensional coordinate under an IMU coordinate system; r is R lr ，R ri Respectively representing rotation matrixes of a right camera to a left camera and a left camera to an IMU coordinate system; t (T) lr ，T ri The translation vectors of the right camera to the left camera and the left camera to the IMU coordinate system, respectively.
Further, the step 3 specifically includes: according to the binocular camera calibrated in the step 2, converting the light center coordinate into a corresponding coordinate under the binocular camera coordinate, and obtaining a light plane equation of the left camera by utilizing the principle of least square fitting plane:
a i x+b i y+c i z+d i ＝0
wherein (x, y, z) is the coordinates of any point in the plane; a, a i 、b i 、c i 、d i To be a constant, a i 、b i 、c i Is the normal vector of the light plane of the left camera and is not 0, d at the same time i Is the intercept of the left camera light plane; i is the number of stripe lasers; the light plane equation of the right camera can be obtained by the same method: a' i x+b′ i y+c′ i z+d′ i ＝0。
Further, the step 4 specifically includes: assuming that the current time stamp is i, since the frequency of the IMU is higher than the update frequency of the point cloud data of the binocular camera, a large amount of IMU data is read in before the next frame of point cloud, namely the time stamp j arrives, the state of the IMU is estimated according to the IMU data in the time interval from i to j, meanwhile, the IMU measurement data in the time interval is utilized to perform pre-integration operation, and the corresponding IMU pre-integration items are as follows:
wherein Δt is the time interval between two consecutive IMU measurements; g w Is a gravity vector in the world coordinate system; r is R k A rotation matrix that turns the IMU coordinate system to the world coordinate system.
Further, the step 5 specifically includes: distortion correction is carried out on input binocular camera image data according to binocular camera internal parameters; preprocessing the corrected left image and the corrected right image, and respectively converting the left image and the corrected right image into a left gray image and a right gray image; the corrected picture is subjected to graying, binarization and noise reduction, and the image convolution operation is carried out on the binarized picture by using Gaussian smoothing operation to remove noise;
any point (x, y) on the image is represented by a Hessian matrix, the matrix consists of four partial derivatives of a certain point of the image along different directions, the matrix can determine the normal direction of the pixel point, and the Taylor expansion fitting sub-pixel coordinates on the normal direction; the eigenvector corresponding to the largest eigenvalue of the Hessian matrix corresponds to the normal direction of the line, using (n x ,n y ) Expressed in terms of points (x 0 ,y 0 ) Setting a threshold for reference points, screening to obtain sub-pixel coordinates of laser center points in left and right images, and constructing a left feature point set (x) l ,y) k Right feature point set (x) r ,y) k ；
The parallax of line laser imaging on the horizontal polar line of the y-th row is calculated by the following formula:
d yk ＝x l -x r ,k＝1,2,...,M
wherein ,dyk Parallax of a kth laser imaging point of an epipolar line with an ordinate y;
calculating to obtain three-dimensional coordinates (X, Y, Z) of each point in space according to an internal reference matrix D and an external reference matrix [ R|T ] of the binocular camera:
wherein f is the focal length of the binocular camera, and the baseline of the binocular camera is given by calibration information; x is x l -x r A disparity between the left view and the right view for a spatial point, given by matching information; (x, y) is the coordinates of the laser spot on the imaging plane.
Further, in the step 6, the transformation matrix is used to update the source point cloud P to obtain a new point cloud P', the new point cloud is used as a new spliced source point cloud, the source point cloud is transformed to the target point cloud, and finally the point cloud splicing is completed:
P′＝RP+T。
further, the step 7 specifically includes:
step 7.1: creating a three-dimensional voxel grid from the input point cloud data, and then in each voxel, using the centers of gravity of all points to approximate other points in the voxel, wherein all points in the voxel use a center of gravity point approximate table;
step 7.2: determining a maximum value max_p and a minimum value min_p of three coordinate axes according to input point cloud data;
step 7.3: determining the side length of a minimum bounding box according to the distance between the maximum value and the minimum value in the range of three coordinate axes, setting the side length of each small grid as L, dividing the number of small grids uniformly divided by the side length of the minimum bounding box of the three coordinate axes into M, N, O, and increasing the required division number of each axis by 1 in order to prevent point cloud data from appearing at the boundary of the bounding box;
step 7.4: for any point P in the point cloud i (X i ,Y i ,Z i ) Calculating index numbers (i, j, k) along the length, width and height, and calculating a grid number h to which the point belongs according to the index numbers; dividing all the point cloud data into small grid cubes, deleting grids without the point cloud data, recalculating the gravity centers of all the points with the same number, replacing other points with the gravity centers, reducing the number of the point clouds, simultaneously keeping the geometric characteristics of the point clouds, and improving the processing speed;
the invention has the beneficial effects that:
the invention uses multi-line structured light to scan the underwater target with high precision, and effectively solves the defects of point cloud distortion, information loss and the like of three-dimensional reconstruction of lattice laser and single-line laser. The underwater target is scanned by the binocular camera with the multi-line laser transmitters and the narrow-band optical filters, two groups of photos are not required to be shot underwater for processing, and scanning efficiency is improved. According to the invention, through the binocular multi-line underwater three-dimensional reconstruction system arranged on the underwater robot, the IMU is utilized to perform pre-integration, the object to be measured is scanned in the moving process, each frame of photo of the underwater binocular camera comprises a plurality of lasers, so that the point cloud density is improved, the scanning frequency is reduced, and the synchronous three-dimensional reconstruction of the underwater environment is realized.
Drawings
FIG. 1 is a schematic representation of the morphology of the process of the present invention.
Fig. 2 is a flow chart of the present invention.
FIG. 3 is a binocular camera epipolar constraint model in the present invention.
Fig. 4 is a flow chart of IMU calibration in the present invention.
Detailed Description
The invention is further described below with reference to the accompanying drawings.
An underwater three-dimensional reconstruction method based on binocular multi-line structured light comprises the steps that two underwater cameras 1 and a multi-line laser transmitter 2 are placed on the same horizontal plane and are fixed in relative positions, and the two underwater cameras 1 and an inertial measurement unit 3 (IMU) are connected with a central processing unit 4.
The platform is arranged at the bottom of the underwater robot, scans towards the sea bottom, and the two underwater cameras 1 incline to the multi-line laser transmitter 2 by a certain angle and are symmetrically arranged, and the two cameras keep a relatively fixed pose relation; the multi-line laser transmitter 2 is used for generating a plurality of 520nm word line structure lights, and the interval distances among the structure lights are required to be consistent; the two underwater cameras 1 are additionally provided with 520nm narrow-band filters, and light signals are allowed to pass through the filters in specific wave bands.
An underwater three-dimensional reconstruction method based on binocular multi-line structured light comprises the following steps:
s1, calibrating parameters of a binocular camera and parameters of an IMU, and calibrating parameters of the binocular camera and the IMU;
s2, calibrating light planes of the left camera and the right camera;
s3, aligning the binocular camera with an IMU time stamp, and collecting synchronous data at the time t; pre-integrating the IMU, processing a laser image of the binocular camera, and extracting point cloud information;
s4, respectively acquiring a left image and a right image of the object to be measured through a left camera and a right camera of the binocular camera, and carrying out three-dimensional correction on the left image and the right image to align the corrected left image and right image rows; matching the corrected left image and right image to obtain a line laser imaging matching point pair; and obtaining left-right image parallax according to the line laser imaging matching point pairs, and calculating three-dimensional data of each laser line according to the left-right image parallax.
S5, performing pose matching correction based on multi-frame registration of the three-dimensional point cloud data, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame registration;
s6, processing the redundant three-dimensional data, and removing noise points and outliers by adopting a voxel filtering method to reconstruct the space three-dimensionally.
The step S1 comprises the following steps:
s11, carrying out three-dimensional calibration on a left camera and a right camera of the binocular camera, carrying out angular point identification calibration on images, namely calibrating a lens focal length f, a base line distance B and an internal reference K of the left camera and the right camera of the binocular camera l 、K r Distortion parameter [ k ] 1 k 2 p 1 p 2 k 3 ]，k 1 k 2 k 3 Radial distortion coefficient, p 1 p 2 Is the tangential distortion coefficient; the external parameters of the binocular camera comprise a rotation matrix R and a translation matrix T;
s12, calibrating the IMU, and acquiring white noise and zero offset of the gyroscope and the accelerometer. Measurement model of accelerometer and gyroscope:
a B ＝T a K a (a S +b a +ν a )
w B ＝T g K g (w S +b g +ν g )
wherein a represents an accelerometer, g represents a gyroscope, B represents an orthogonal reference coordinate system, S represents a non-orthogonal tuning coordinate, T represents a transformation matrix of an axis deviation, K represents a scale factor, a S 、w S Representing true values, b, v representing bias and white noise.
Issuing an IMU topic through an ROS command, standing the IMU still to record an IMU data packet, calibrating a program to play back the recorded data packet, and calculating white noise, bias and average value of the triaxial of the accelerometer and the gyroscope;
s13, calibrating the left and right cameras and the IMU external parameters of the binocular camera means that the IMU coordinate system is set as a world coordinate system, and then the conversion relation from the image points of the left and right cameras of the binocular camera to the IMU coordinate system is as follows:
wherein And->Is two-dimensional coordinates under the coordinates of the left camera and the right camera, < >>Is three-dimensional coordinates in an IMU coordinate system, R lr ，R ri Respectively representing rotation matrixes of a right camera to a left camera and a left camera to an IMU coordinate system; t (T) lr ，T ri Translation vectors from the right camera to the left camera and from the left camera to the IMU coordinate system;
the step S2 comprises the following steps:
s21, calibrating the light planes of the left and right cameras, calibrating the external parameters of the underwater camera in the step S2, converting the light center coordinates into corresponding coordinates under a camera coordinate system, and obtaining a light plane equation of the left camera by utilizing the principle of least square fitting of the plane:
a i x+b i y+c i z+d i ＝0
wherein x, y, z are the coordinates of any point in the plane, a i 、b i 、c i 、d i To be a constant, a i 、b i 、c i Is the normal vector of the light plane of the left camera and is not 0, d at the same time i The intercept of the light plane of the left camera is i, and the number of stripe lasers is i; the right camera light plane equation can be obtained by the same method: a' i x+b′ i y+c′ i z+d′ i ＝0。
The binocular vision and IMU-based three-dimensional reconstruction method of an underwater scene according to claim 2, wherein:
the step S3 comprises the following steps:
s31, assuming that the current timestamp is i, since the frequency of the IMU is generally far higher than the update frequency of the point cloud data of the binocular camera, a large amount of IMU data can be read in before the next frame of point cloud, namely the timestamp j, is arrived, the state of the IMU can be estimated according to the IMU data in the time interval from i to j, and meanwhile, the IMU measurement data in the time interval is utilized for carrying out pre-integration operation;
the corresponding IMU pre-integral term is:
wherein: Δt is the time interval between two consecutive IMU measurements; g w Is a gravity vector in the world coordinate system; r is R k A rotation matrix for converting the IMU coordinate system to a world coordinate system;
the step S4 comprises the following steps:
s41, carrying out distortion correction on input binocular camera image data according to the binocular camera internal parameters; preprocessing the corrected left image and the corrected right image, and respectively converting the left image and the corrected right image into a left gray image and a right gray image; the corrected picture is subjected to graying, binarization and noise reduction, and the image convolution operation is carried out on the binarized picture by using Gaussian smoothing operation to remove noise;
s42, any point (x, y) on the image can be represented by a Hessian matrix, the matrix consists of four partial derivatives of a certain point of the image along different directions, and the matrix can determine the normal direction of the pixel point; the normal taylor expansion fits the sub-pixel coordinates. The eigenvector corresponding to the largest eigenvalue of the Hessian matrix corresponds to the normal direction of the line, using (n x ,n y ) Expressed in terms of points (x 0 ,y 0 ) Setting a threshold for reference points, screening to obtain sub-pixel coordinates of laser center points in left and right images, and constructing a left feature point set (x) l ,y) k Right feature point set (x) r ,y) k 。
S43, calculating parallax of line laser imaging on a y-th horizontal polar line according to the following formula:
d yk ＝x l -x r ,k＝1,2,...,M
wherein ,dyk Parallax of a kth laser imaging point of an epipolar line with an ordinate y;
s44, calculating to obtain three-dimensional coordinates (X, Y, Z) of each point in space according to the internal reference matrix D and the external reference matrix [ R|T ] of the binocular camera and the following formula
Wherein f is an extrinsic matrix [ R|T ]]The focal length of the binocular camera, the base line of the binocular camera, is given by calibration information; x is x l -x r A disparity between the left view and the right view for a spatial point, given by matching information; (x, y) is the coordinates of the laser spot on the imaging plane.
The step S5 comprises the following steps:
s51, updating the source point cloud P by using the transformation matrix to obtain a new point cloud P', taking the new point cloud as a new spliced source point cloud, transforming the source point cloud to a target point cloud, and finally finishing point cloud splicing.
P′＝RP+T
The step S6 comprises the following steps:
s61, creating a three-dimensional voxel grid from the input point cloud data, and then in each voxel, using the centers of gravity of all points to approximate and display other points in the voxel, wherein all points in the voxel use a center of gravity point approximate table;
s62, determining a maximum value max_p and a minimum value min_p of three coordinate axes according to input point cloud data;
s63, determining the side length of a minimum bounding box according to the distance between the maximum value and the minimum value in the range of three coordinate axes, setting the side length of each small grid as L, dividing the number of small grids uniformly divided by the side length of the minimum bounding box of the three coordinate axes into M, N, O, and increasing the required division number of each axis by 1 in order to prevent point cloud data from appearing at the bounding box boundary;
s64, for any point P in the point cloud i (X i ,Y i ,Z i ) Calculating index numbers (i, j, k) along the length, width and height, and calculating a grid number h to which the point belongs according to the index numbers; dividing all the point cloud data into small grid cubes, deleting grids without the point cloud data, recalculating the gravity centers of all the points with the same number, replacing other points with the gravity centers, reducing the number of the point clouds, simultaneously keeping the geometric characteristics of the point clouds, and improving the processing speed.
h＝i+j×M+k×M×N
Example 1:
the invention discloses a binocular multi-line structured light underwater three-dimensional reconstruction system, which comprises two cameras, a 520nm multi-line laser transmitter, a central processing unit and an Inertial Measurement Unit (IMU), wherein the two cameras are additionally provided with 520nm narrow-band filters, and the relative positions of the underwater cameras and a laser light source are fixed; firstly, calibrating a binocular camera, an IMU internal parameter and an IMU external parameter; secondly, aligning the binocular camera with an IMU timestamp, carrying out synchronous data acquisition at the time t, pre-integrating the IMU, preprocessing an image of the binocular camera, extracting laser centers to construct left and right characteristic points, obtaining matching point pairs in the left and right characteristic point sets by using polar constraint, and calculating three-dimensional data of each laser line according to left and right image parallax; performing pose matching correction based on multi-frame matching, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame matching; and processing the redundant three-dimensional data, removing noise points and outliers by adopting a voxel filtering method, realizing synchronous three-dimensional reconstruction of the underwater environment, improving the density of point cloud and reducing the scanning frequency.
The invention comprises two underwater cameras 1 added with 520nm narrow-band filters, a 520nm multi-line laser transmitter 2, an inertial measurement unit 3 (IMU) and a central processing unit 4. The binocular underwater camera and the multi-line laser transmitter are fixed in relative positions, and the central processing unit is connected with the binocular camera in a USB mode.
The whole three-dimensional scanning system is arranged at the bottom of the underwater robot and scans along with the movement of the underwater robot;
an underwater three-dimensional reconstruction method based on binocular multi-line structured light is carried out according to the following steps:
s1, an underwater scanning system is arranged at the bottom of an underwater robot, scans towards the sea bottom, and two underwater cameras 1 are inclined at a certain angle towards a multi-line laser transmitter 2 and are symmetrically arranged, and a relatively fixed pose relation is kept between the two cameras; the multi-line laser transmitter 2 is used for generating a plurality of 520nm word line structure lights, and the interval distances among the structure lights are required to be consistent; the underwater camera 1 is additionally provided with a 520nm narrow-band filter, and allows light signals to pass through in a specific wave band.
S2, calibrating internal parameters of the underwater binocular camera, IMU parameters and external parameters of the binocular camera and the IMU;
s3, calibrating light planes of the left camera and the right camera;
s4, aligning the binocular camera with an IMU time stamp, and carrying out synchronous data acquisition at the time t to pre-integrate the IMU; the whole scanning system moves along with the advancing of the underwater robot;
s5, the central processing unit respectively acquires a left image and a right image of the object to be detected through a left camera and a right camera of the binocular camera, and performs three-dimensional correction on the left image and the right image to align the corrected left image and right image; matching the corrected left image and right image to obtain a line laser imaging matching point pair; and obtaining left-right image parallax according to the line laser imaging matching point pairs, and calculating three-dimensional data of each laser line according to the left-right image parallax.
S6, performing pose matching correction based on multi-frame matching of the three-dimensional point cloud data, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame matching;
and S7, processing the redundant three-dimensional data, and removing noise points and outliers by adopting a voxel filtering method to reconstruct the space three-dimensionally.
S2 comprises the following steps:
performing three-dimensional calibration on a left camera and a right camera of the binocular camera, and performing corner recognition calibration on images to obtain lens focal lengths, base line distances, internal parameters and distortion parameters of the left and right cameras of the binocular camera; wherein the left camera is internally provided with reference K l Internal reference K of right camera r The method comprises the following steps of:
wherein ,f xl f yl f xr f yr scale factors on the x and y axes, also called normalized focal length; (c) xl ,c yl )、(c xr ,c yr ) Is the optical center coordinates of the left and right cameras.
The external parameters of the binocular camera, namely the rotation matrix R, and the translation vector T are:
wherein rs (s=1, 2., 9) is the first of the rotation matrices Rs elements, tx, ty and tz are three components of the translation vector T, respectively;
distortion parameter [ k ] 1 ,k 2 ,p 1 ,p 2 ,k 3 ]，k 1 、k 2 、k 3 Radial distortion coefficient, p 1 、p 2 Is the tangential distortion coefficient;
IMU calibration as shown in fig. 3, the initialization time T during calibration. The time t for which the rotation is kept still is generally 36 to 50 times. In general, the more and better this number is, the fewer rotations must be, the greater the number of parameters required to solve, since this avoids singularities in the solution.
External parameter calibration of a binocular camera and an IMU:
the camera is used to capture images containing the april grid calibration plate while the camera is fully moved along the 6 axes of the IMU, including straight line movement and rotation along the accelerometer axis, and a free movement is performed to cover IMU information in all directions.
The conversion relation from the image points of the left camera and the right camera to the IMU coordinate system of the platform is as follows:
wherein And->Is two-dimensional coordinates under the coordinates of the left camera and the right camera, < >>Is three-dimensional coordinates in an IMU coordinate system, R lr ，R ri Respectively representing rotation matrixes of a right camera to a left camera and a left camera to an IMU coordinate system; t (T) lr ，T ri Translation vectors from the right camera to the left camera and from the left camera to the IMU coordinate system;
calibrating the light planes of the left and right cameras, calibrating the external parameters of the underwater camera in step S2, converting the light center coordinates into corresponding coordinates under the camera coordinate system, and obtaining a light plane equation of the left camera by utilizing the principle of least square fitting of the plane: a, a i x+b i y+c i z+d i =0, where x, y, z are the coordinates of any point in this plane; a, a i 、b i 、c i 、d i To be a constant, a i 、b i 、c i Is the normal vector of the light plane of the left camera and is not 0, d at the same time i The intercept of the light plane of the left camera is i, and the number of stripe lasers is i; the right camera light plane equation can be obtained by the same method: a, a i ′x+b i ′y+c i ′z+d i ′＝0。
The central processing unit controls the binocular camera to shoot pictures with a plurality of lasers, unifies the time stamp of the binocular camera and the IMU, pre-integrates the IMU, extracts laser point cloud information from the binocular image, reads the data of the IMU, and acquires an inter-frame rotation matrix and a translation vector.
The state information of the AUV can be further solved by the acceleration and the angular velocity, and the state vector at a certain moment is defined as follows:
wherein ,ba ,b g The bias of the accelerometer and gyroscope, respectively, are considered herein to be constant,the position and the speed of the machine body coordinate system relative to the inertial coordinate system at the moment k can be obtained by integrating acceleration values>The orientation of the underwater robot coordinate system is defined, expressed in terms of quaternions, which are related to angular velocity. The quaternion q is derived from time to obtain:
the integral of angular velocity is the rotation angle, the integral of acceleration is the velocity, and the integral of velocity is the displacement. The relative pose change during this time can be known by IMU integration over two frames of time t. And further, the position of each moment can be obtained in an accumulated manner.
wherein ,ΔRij Representing a relative rotation between the i-th image and the j-th image;a measurement representing the rotation of the kth gyroscope; />Indicating zero offset corresponding to rotation measurement of the kth gyroscope; Δt represents the time interval between two gyroscope measurements; deltaV ij Representing a change in the relative degree between the ith image and the jth image; />Representing the measurement of the kth accelerometer; />Indicating zero offset corresponding to the kth accelerometer measurement; ΔP ij Representing the relative position change between the i-th image and the j-th image.
The central processing unit carries out distortion correction processing on the photos of the left and right cameras, and as the camera with the 520nm laser emitter and the 520nm narrow-band filter is adopted, the background noise of the underwater photo is reduced, the collected photo can be directly subjected to graying, binarization and noise reduction, and a left and right characteristic point set is constructed by extracting the left and right image laser centers by using a steger algorithm.
Graying: the process of converting a color image into a grayscale image becomes a graying process of the image. The calculation amount of the subsequent image is reduced, and the gray image can still reflect the distribution and characteristics of the chromaticity and brightness level of the whole and part of the photo:
Gray＝0.299R+0.587G+0.114B
histogram equalization: the method is mainly used for enhancing the contrast of the image with smaller dynamic range, the histogram of the original image is transformed into a uniformly distributed form, and the dynamic range of the gray value of the pixel is enhanced, so that the effect of enhancing the overall contrast of the image is achieved. The gray level histogram is a method for counting the occurrence frequency of all pixels in the digital image according to the gray level value:
where k is the k-th gray scale value of the image f (m, n), n k The number of pixels with a gray value k in f (m, n), n being the total number of pixels of the image, and L being the number of gray levels.
Adaptive threshold image binarization: according to the brightness distribution of different areas of the image, the local threshold value is calculated, and for the different areas of the image, different threshold values can be calculated adaptively.
Gaussian smoothing filtering: and carrying out weighted average on the whole image, wherein the value of each pixel point is obtained by carrying out weighted average on the pixel point and other pixel values in the neighborhood. Each pixel in the image is scanned by a convolution, and the value of the convolutionally centered pixel point is replaced by the convolutionally determined weighted average gray value of the pixels in the neighborhood.
Any point (x, y) on the image, the Hessian matrix can be expressed as:
the matrix consists of four partial derivatives of a certain point of an image along different directions, and the normal direction of the pixel point can be determined by the matrix; the normal taylor expansion fits the sub-pixel coordinates. The eigenvector corresponding to the largest eigenvalue of the Hessian matrix corresponds to the normal direction of the line, using (n x ,n y ) Expressed in terms of points (x 0 ,y 0 ) Is a datum point; setting a threshold value for screening, wherein the coordinates of the sub-pixels at the center of the line are as follows:
obtaining coordinates of laser center points in the left and right images, and constructing a left feature point set (x l ,y) k Right feature point set (x) r ,y) k 。
The parallax of line laser imaging on the horizontal polar line of the y-th row is calculated by the following formula:
d yk ＝x l -x r ,k＝1,2,...,M
wherein ,dyk Parallax of a kth laser imaging point of an epipolar line with an ordinate y;
and calculating to obtain three-dimensional coordinates (X, Y, Z) of each point in space according to an internal reference matrix D and an external reference matrix [ R|T ] of the binocular camera and the following formula:
wherein f is an extrinsic matrix [ R|T ]]The focal length of the medium binocular camera, B, is the distance between the left camera and the right camera in the internal reference matrix D, and is given by calibration information; x is x l -x r A disparity between the left view and the right view for a spatial point, given by matching information; (x, y) is the coordinates of the laser spot on the imaging plane.
Performing pose matching correction based on multi-frame matching of three-dimensional point cloud data, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame matching:
and updating the source point cloud P by using the transformation matrix obtained by solving to obtain a new point cloud P' =RP+T, taking the new point cloud as a new spliced source point cloud, transforming the source point cloud to a target point cloud, and finally finishing point cloud splicing.
Processing redundant three-dimensional data, removing noise points and outliers by adopting a voxel filtering method, and performing space three-dimensional reconstruction
Creating a three-dimensional voxel grid from the input point cloud data, and then in each voxel, using the centers of gravity of all points to approximate other points in the voxel, wherein all points in the voxel use a center of gravity point approximate table;
determining a maximum value max_p and a minimum value min_p of three coordinate axes according to input point cloud data; determining the side length of a minimum bounding box according to the distance between the maximum value and the minimum value in the range of three coordinate axes, setting the side length of each small grid as L, dividing the number of small grids uniformly divided by the side length of the minimum bounding box of the three coordinate axes into M, N, O, and increasing the required division number of each axis by 1 in order to prevent point cloud data from appearing at the boundary of the bounding box;
for any point P in the point cloud i (X i ,Y i ,Z i ) Calculating index numbers (i, j, k) along the length, width and height, and calculating a grid number h to which the point belongs according to the index numbers; dividing all the point cloud data into small grid cubes, deleting grids without the point cloud data, recalculating the gravity centers of all the points with the same number, replacing other points with the gravity centers, reducing the number of the point clouds, simultaneously keeping the geometric characteristics of the point clouds, and improving the processing speed.
h＝i+j×M+k×M×N
wherein ,
the invention has the beneficial effects that:
the invention uses multi-line structured light to scan the underwater target with high precision, and effectively solves the defects of point cloud distortion, information loss and the like of three-dimensional reconstruction of lattice laser and single-line laser. The underwater target is scanned by the binocular camera with the multi-line laser transmitters and the narrow-band optical filters, two groups of photos (one group is provided with laser and the other group is not provided with laser) are not required to be shot underwater for processing, and the scanning efficiency is improved. The binocular multi-line underwater three-dimensional reconstruction system installed on the underwater robot scans the object to be detected in the moving process by pre-integrating the IMU, and each frame of photo of the underwater camera contains a plurality of lasers, so that the point cloud density is improved, the scanning frequency is reduced, and the synchronous three-dimensional reconstruction of the underwater environment is realized.
The above description is only of the preferred embodiments of the present invention and is not intended to limit the present invention, but various modifications and variations can be made to the present invention by those skilled in the art. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention should be included in the protection scope of the present invention.
Claims (7)
1. An underwater three-dimensional reconstruction method based on binocular multi-line structured light is characterized by comprising the following steps:
step 1: a binocular camera and a multi-line laser transmitter are arranged at the bottom of the underwater robot; the binocular camera comprises a left camera and a right camera which are respectively arranged at two sides of the multi-line laser transmitter;
step 2: calibrating the binocular camera and an inertial measurement unit IMU;
step 3: calibrating light planes of the left camera and the right camera;
step 4: aligning the binocular camera with an inertial measurement unit IMU time stamp, and carrying out synchronous data acquisition at the time t; pre-integrating an Inertial Measurement Unit (IMU), processing a laser image of a binocular camera, and extracting point cloud information;
step 5: the method comprises the steps that a left camera and a right camera are used for respectively obtaining a left image and a right image of an object to be measured, and the left image and the right image are subjected to three-dimensional correction, so that corrected left image and corrected right image are aligned; matching the corrected left image and right image to obtain a line laser imaging matching point pair; obtaining left and right image parallax according to the line laser imaging matching point pairs, and calculating three-dimensional data of each laser line according to the left and right image parallax;
step 6: performing pose matching correction based on multi-frame registration of the three-dimensional point cloud data, and obtaining multi-frame three-dimensional point cloud data by using a rotation matrix R and a translation vector T after inter-frame registration;
step 7: and processing the redundant three-dimensional data, and removing noise points and outliers by adopting a voxel filtering method to reconstruct the space three-dimensionally.
2. The underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: the step 2 of calibrating the binocular camera comprises calibrating the lens focal length f of the binocular camera, the base line distance B and the internal reference K of the left camera l Internal reference K of right camera r Distortion parameter [ k ] 1 k 2 p 1 p 2 k 3 ]A rotation matrix R and a translation matrix T of the binocular camera; wherein k is 1 、k 2 、k 3 For radial distortion coefficient, p 1 、p 2 Is a tangential distortion coefficient;
the calibrating the inertial measurement unit IMU comprises setting an IMU coordinate system as a world coordinate system, and then converting the binocular camera coordinate system into the IMU coordinate system into the following relationship:
wherein ,and->Is two-dimensional coordinates under binocular camera coordinates, +.>Is a three-dimensional coordinate under an IMU coordinate system; r is R lr ，R ri Respectively representing rotation matrixes of a right camera to a left camera and a left camera to an IMU coordinate system; t (T) lr ，T ri The translation vectors of the right camera to the left camera and the left camera to the IMU coordinate system, respectively.
3. The underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: the step 3 specifically comprises the following steps: according to the binocular camera calibrated in the step 2, converting the light center coordinate into a corresponding coordinate under the binocular camera coordinate, and obtaining a light plane equation of the left camera by utilizing the principle of least square fitting plane:
a i x+b i y+c i z+d i ＝0
wherein (x, y, z) is the coordinates of any point in the plane; a, a i 、b i 、c i 、d i To be a constant, a i 、b i 、c i Is the normal vector of the light plane of the left camera and is not 0, d at the same time i Is the intercept of the left camera light plane; i is the number of stripe lasers; the light plane equation of the right camera can be obtained by the same method: a' i x+b′ i y+c′ i z+d′ i ＝0。
4. The underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: the step 4 specifically comprises the following steps: assuming that the current time stamp is i, since the frequency of the IMU is higher than the update frequency of the point cloud data of the binocular camera, a large amount of IMU data is read in before the next frame of point cloud, namely the time stamp j arrives, the state of the IMU is estimated according to the IMU data in the time interval from i to j, meanwhile, the IMU measurement data in the time interval is utilized to perform pre-integration operation, and the corresponding IMU pre-integration items are as follows:
wherein Δt is the time interval between two consecutive IMU measurements; g w Is a gravity vector in the world coordinate system; r is R k A rotation matrix that turns the IMU coordinate system to the world coordinate system.
5. The underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: the step 5 specifically comprises the following steps: distortion correction is carried out on input binocular camera image data according to binocular camera internal parameters; preprocessing the corrected left image and the corrected right image, and respectively converting the left image and the corrected right image into a left gray image and a right gray image; the corrected picture is subjected to graying, binarization and noise reduction, and the image convolution operation is carried out on the binarized picture by using Gaussian smoothing operation to remove noise;
any point (x, y) on the image is represented by a Hessian matrix, the matrix consists of four partial derivatives of a certain point of the image along different directions, the matrix can determine the normal direction of the pixel point, and the Taylor expansion fitting sub-pixel coordinates on the normal direction; the eigenvector corresponding to the largest eigenvalue of the Hessian matrix corresponds to the normal direction of the line, using (n x ,n y ) Expressed in terms of points (x 0 ,y 0 ) Setting a threshold for reference points, screening to obtain sub-pixel coordinates of laser center points in left and right images, and constructing a left feature point set (x) l ,y) k Right feature point set (x) r ,y) k ；
The parallax of line laser imaging on the horizontal polar line of the y-th row is calculated by the following formula:
d yk ＝x l -x r ,k＝1,2,...,M
wherein ,dyk Parallax of a kth laser imaging point of an epipolar line with an ordinate y;
calculating to obtain three-dimensional coordinates (X, Y, Z) of each point in space according to an internal reference matrix D and an external reference matrix [ R|T ] of the binocular camera:
wherein f is the focal length of the binocular camera, and the baseline of the binocular camera is given by calibration information; x is x l -x r A disparity between the left view and the right view for a spatial point, given by matching information; (x, y) is the coordinates of the laser spot on the imaging plane.
6. The underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: in the step 6, the transformation matrix is used to update the source point cloud P to obtain a new point cloud P', the new point cloud is used as a new spliced source point cloud, the source point cloud is transformed to the target point cloud, and finally the point cloud splicing is completed:
P′＝RP+T。
7. the underwater three-dimensional reconstruction method based on binocular multi-line structured light according to claim 1, wherein the method comprises the following steps: the step 7 specifically comprises the following steps:
step 7.1: creating a three-dimensional voxel grid from the input point cloud data, and then in each voxel, using the centers of gravity of all points to approximate other points in the voxel, wherein all points in the voxel use a center of gravity point approximate table;
step 7.2: determining a maximum value max_p and a minimum value min_p of three coordinate axes according to input point cloud data;
step 7.3: determining the side length of a minimum bounding box according to the distance between the maximum value and the minimum value in the range of three coordinate axes, setting the side length of each small grid as L, dividing the number of small grids uniformly divided by the side length of the minimum bounding box of the three coordinate axes into M, N, O, and increasing the required division number of each axis by 1 in order to prevent point cloud data from appearing at the boundary of the bounding box;
step 7.4: for any point P in the point cloud i (X i ,Y i ,Z i ) Calculating index numbers (i, j, k) along the length, width and height, and calculating a grid number h to which the point belongs according to the index numbers; all points are takenDividing the cloud data into small grid cubes, deleting grids without point cloud data, recalculating the gravity centers of all points in the same number, replacing other points with the gravity center points, reducing the number of the point cloud, simultaneously keeping the geometric characteristics of the point cloud, and improving the processing speed;
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310948944.0A CN116977391A (en) | 2023-07-31 | 2023-07-31 | Underwater three-dimensional reconstruction method based on binocular multi-line structured light |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310948944.0A CN116977391A (en) | 2023-07-31 | 2023-07-31 | Underwater three-dimensional reconstruction method based on binocular multi-line structured light |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN116977391A true CN116977391A (en) | 2023-10-31 |
Family
ID=88479132
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310948944.0A Pending CN116977391A (en) | 2023-07-31 | 2023-07-31 | Underwater three-dimensional reconstruction method based on binocular multi-line structured light |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116977391A (en) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117173342A (en) * | 2023-11-02 | 2023-12-05 | 中国海洋大学 | Underwater monocular and binocular camera-based natural light moving three-dimensional reconstruction device and method |
| CN117994446A (en) * | 2024-04-07 | 2024-05-07 | 华东交通大学 | Light fusion complementary three-dimensional reconstruction method and system based on polarized binocular line structure |
-
2023
- 2023-07-31 CN CN202310948944.0A patent/CN116977391A/en active Pending
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117173342A (en) * | 2023-11-02 | 2023-12-05 | 中国海洋大学 | Underwater monocular and binocular camera-based natural light moving three-dimensional reconstruction device and method |
| CN117994446A (en) * | 2024-04-07 | 2024-05-07 | 华东交通大学 | Light fusion complementary three-dimensional reconstruction method and system based on polarized binocular line structure |
| CN117994446B (en) * | 2024-04-07 | 2024-06-11 | 华东交通大学 | Light fusion complementary three-dimensional reconstruction method and system based on polarized binocular line structure |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109615652B (en) | Depth information acquisition method and device | |
| CN110264567B (en) | Real-time three-dimensional modeling method based on mark points | |
| CN116977391A (en) | Underwater three-dimensional reconstruction method based on binocular multi-line structured light | |
| CN109166077B (en) | Image alignment method and device, readable storage medium and computer equipment | |
| Sedlazeck et al. | 3d reconstruction based on underwater video from rov kiel 6000 considering underwater imaging conditions | |
| CN114004941A (en) | Indoor scene three-dimensional reconstruction system and method based on nerve radiation field | |
| US5259037A (en) | Automated video imagery database generation using photogrammetry | |
| CN111369630A (en) | Method for calibrating multi-line laser radar and camera | |
| CN109544679A (en) | The three-dimensional rebuilding method of inner wall of the pipe | |
| CN110334701B (en) | Data acquisition method based on deep learning and multi-vision in digital twin environment | |
| CN109712232B (en) | Object surface contour three-dimensional imaging method based on light field | |
| CN107170037A (en) | A kind of real-time three-dimensional point cloud method for reconstructing and system based on multiple-camera | |
| CN205451195U (en) | Real -time three -dimensional some cloud system that rebuilds based on many cameras | |
| CN108663026A (en) | A kind of vibration measurement method | |
| CN112561996A (en) | Target detection method in autonomous underwater robot recovery docking | |
| CN112929626A (en) | Three-dimensional information extraction method based on smartphone image | |
| CN114359406A (en) | Calibration of auto-focusing binocular camera, 3D vision and depth point cloud calculation method | |
| CN113034571A (en) | Object three-dimensional size measuring method based on vision-inertia | |
| CN110021035B (en) | Marker of Kinect depth camera and virtual marker tracking method based on marker | |
| CN112465778A (en) | Underwater fish shoal observation device and method | |
| CN116205961A (en) | Automatic registration method and system for multi-lens combined image and laser radar point cloud | |
| CN113345084B (en) | Three-dimensional modeling system and three-dimensional modeling method | |
| CN112712566B (en) | Binocular stereo vision sensor measuring method based on structure parameter online correction | |
| CN117406234A (en) | Target ranging and tracking method based on single-line laser radar and vision fusion | |
| CN116912393A (en) | Face reconstruction method and device, electronic equipment and readable storage medium |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |