WO2017216999A1 - 波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム - Google Patents

波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム Download PDF

Info

Publication number
WO2017216999A1
WO2017216999A1 PCT/JP2017/002787 JP2017002787W WO2017216999A1 WO 2017216999 A1 WO2017216999 A1 WO 2017216999A1 JP 2017002787 W JP2017002787 W JP 2017002787W WO 2017216999 A1 WO2017216999 A1 WO 2017216999A1
Authority
WO
WIPO (PCT)
Prior art keywords
wave source
function
cross
input signal
calculating
Prior art date
Application number
PCT/JP2017/002787
Other languages
English (en)
French (fr)
Inventor
正徳 加藤
裕三 仙田
Original Assignee
日本電気株式会社
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by 日本電気株式会社 filed Critical 日本電気株式会社
Priority to JP2018523292A priority Critical patent/JP6874766B2/ja
Priority to US16/308,860 priority patent/US11454694B2/en
Publication of WO2017216999A1 publication Critical patent/WO2017216999A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/802—Systems for determining direction or deviation from predetermined direction
    • G01S3/808—Systems for determining direction or deviation from predetermined direction using transducers spaced apart and measuring phase or time difference between signals therefrom, i.e. path-difference systems
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/802—Systems for determining direction or deviation from predetermined direction
    • G01S3/805—Systems for determining direction or deviation from predetermined direction using adjustment of real or effective orientation of directivity characteristics of a transducer or transducer system to give a desired condition of signal derived from that transducer or transducer system, e.g. to give a maximum or minimum signal
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/80—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves
    • G01S3/86—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using ultrasonic, sonic or infrasonic waves with means for eliminating undesired waves, e.g. disturbing noises

Definitions

  • the present invention relates to a wave source direction estimation device, a wave source direction estimation system, a wave source direction estimation method, and a wave source direction estimation program.
  • Non-Patent Document 1 discloses a technique for estimating the direction of a sound source (sound wave generation source, generation location) after obtaining a plurality of cross-correlation functions using sound reception signals of two or more microphones.
  • a sound source sound wave generation source, generation location
  • cross-correlation functions between received signals are obtained, and the cross-correlation functions corresponding to each combination are summed to calculate one cross-correlation function.
  • the arrival direction of the sound wave is estimated by calculating the time difference at which the cross-correlation function obtained by the sum gives the maximum value as the arrival time difference of the sound wave.
  • the direction of the wave source cannot be estimated with high accuracy in an environment where the ambient noise level is high.
  • An object of the present invention is to provide a technique for solving the above-described problems.
  • a wave source direction estimation device In order to achieve the above object, a wave source direction estimation device according to the present invention is provided.
  • a plurality of input signal acquisition means for acquiring signals generated at the wave source as input signals;
  • Cross-correlation function calculating means for calculating a cross-correlation function based on the input signal acquired by the input signal acquiring means;
  • An envelope function extracting means for extracting an envelope function based on the calculated cross-correlation function;
  • An integrated envelope function calculating means for calculating an integrated envelope function by integrating the extracted envelope functions;
  • estimated direction information generating means Based on the calculated integrated envelope function, estimated direction information generating means for generating estimated direction information of the wave source; Equipped with.
  • a wave source direction estimation system includes: A plurality of input signal acquisition means for acquiring signals generated at the wave source as input signals; Cross-correlation function calculation means for calculating a cross-correlation function based on the input signals acquired by the plurality of input signal acquisition means; An envelope function extracting means for extracting an envelope function based on the calculated cross-correlation function; An integrated envelope function calculating means for calculating an integrated envelope function by integrating the extracted envelope functions; Based on the calculated integrated envelope function, estimated direction information generating means for generating estimated direction information of the wave source; Display means for displaying the estimated direction information; Equipped with.
  • a wave source direction estimation method includes: An input signal acquisition step of acquiring a signal generated at the wave source as an input signal; A cross-correlation function calculating step for calculating a cross-correlation function based on the input signal acquired in the input signal acquiring step; An envelope function extracting step for extracting an envelope function based on the calculated cross-correlation function; An integrated envelope function calculating step of calculating an integrated envelope function by integrating the extracted envelope functions; An estimated direction information generation step for generating estimated direction information of the wave source based on the calculated integrated envelope function; including.
  • a wave source direction estimation program includes: An input signal acquisition step of acquiring a signal generated at the wave source as an input signal; A cross-correlation function calculating step for calculating a cross-correlation function based on the input signal acquired in the input signal acquiring step; An envelope function extracting step for extracting an envelope function based on the calculated cross-correlation function; An integrated envelope function calculating step of calculating an integrated envelope function by integrating the extracted envelope functions; An estimated direction information generation step for generating estimated direction information of the wave source based on the calculated integrated envelope function; Is executed on the computer.
  • the direction of the wave source can be estimated with high accuracy even in an environment where the ambient noise level is high.
  • the “voice signal” is a direct electrical change that occurs in accordance with voice and other sounds, and is used to transmit voice and other sounds, and is not limited to voice.
  • the estimation target of the wave source direction estimation device is not limited to the generation source of sound waves that are vibration waves of air or water. It can also be applied to sources of vibration waves using earth and solids such as earthquakes and landslides. In that case, a vibration sensor is used instead of a microphone as a device that converts vibration waves into electrical signals.
  • the wave source direction estimation apparatus can be applied to the case where the direction is estimated using radio waves as well as gas / liquid / solid vibration waves. In that case, an antenna is used as a device that converts radio waves into electrical signals.
  • the wave source is assumed to be a sound source.
  • a wave source direction estimation apparatus 100 as a first embodiment of the present invention will be described with reference to FIG.
  • the wave source direction estimation device 100 is a device that estimates the direction of the wave source position based on the acquired signal.
  • the wave source direction estimation apparatus 100 includes an input signal acquisition unit 101, a cross-correlation function calculation unit 102, an envelope function extraction unit 103, an integrated envelope function calculation unit 104, and an estimated direction information generation unit 105. Including.
  • the plurality of input signal acquisition units 101 acquire signals generated at the wave source as input signals.
  • the cross-correlation function calculation unit 102 calculates a cross-correlation function based on the input signals acquired by the plurality of input signal acquisition units 101.
  • the envelope function extraction unit 103 extracts an envelope function based on the calculated cross-correlation function.
  • the integrated envelope function calculation unit 104 integrates the extracted envelope functions to calculate an integrated envelope function.
  • the estimated direction information generation unit 105 generates estimated direction information of the wave source based on the calculated integrated envelope function.
  • the direction of the wave source can be estimated with high accuracy even in an environment where the ambient noise level is high.
  • Non-Patent Document 1 and Non-Patent Document 2 it is difficult to accurately estimate the direction of a sound source that exists far away in an environment with a high ambient noise level such as outdoors. For example, when the sound source to be estimated (target sound source) is present at a location far away from the microphone, the volume of sound radiated from the target sound source is significantly reduced when reaching the microphone. For this reason, the sound of the target sound source is buried in the ambient noise, and the direction estimation accuracy of the target sound source may be reduced.
  • target sound source target sound source
  • FIG. 2 is a block diagram for explaining the configuration of the wave source direction estimation apparatus 200 according to the present embodiment.
  • the wave source direction estimation apparatus 200 of this embodiment functions as a part of an apparatus such as a digital video camera, a smartphone, a mobile phone, a notebook computer, or a passive sonar. It is also installed in an abnormal sound detection device that detects abnormalities based on voices and sounds such as suspicious drone detection, scream detection, and vehicle accident detection.
  • application examples of the wave source direction estimation apparatus 200 according to the present embodiment are not limited to these, and can be applied to any wave source direction estimation apparatus that is required to estimate the direction of a target sound source from received sound.
  • the wave source direction estimation apparatus 200 includes a cross-correlation function calculation unit 201, an envelope function extraction unit 202, an integrated envelope function calculation unit 203, an estimated direction information generation unit 204, and a relative delay time calculation unit 206.
  • a sound of a target sound source and a sound in which various noises generated around a microphone (hereinafter referred to as a microphone) that is a sound collecting device are mixed are input as digital signals (sample value series).
  • the sound signal input to the input terminals 20 1 to 20 M is referred to as an input signal in this embodiment.
  • the sound input to the input terminal is collected by a microphone that is a sound collector. Since there are a plurality of input terminals, when collecting the sound of the target sound source, M microphones as many as the number of terminals are used simultaneously. In this embodiment, it is assumed that the input terminal and the microphone correspond one-to-one, and the sound collected by the a-th microphone is supplied to the a-th input terminal. Therefore, the input signal supplied to the a-th input terminal is also referred to as “a-th microphone input signal”.
  • the M microphones that collect the input sound to the wave source direction estimation apparatus 200 are arranged in various layouts.
  • the most basic arrangement is an arrangement arranged in a line on a straight line.
  • lattice form it may arrange
  • it may be arranged on a vertex or a sphere of a regular polyhedron. It is assumed that the input signal supplied to the input terminal is a sound collected by the microphone arranged in this way.
  • the wave source direction estimation apparatus 200 estimates the direction of the sound source using the time difference at which the sound of the target sound source reaches a plurality of microphones. For this reason, since the microphone interval is also important information, not only the input signal but also the microphone position information is supplied to the wave source direction estimation apparatus 200.
  • M types of input signals xa (t) are input to the cross-correlation function calculation unit 201.
  • a is an input terminal number.
  • the cross-correlation function calculation unit 201 cuts out a waveform having an appropriate length from the input signal supplied from the input terminal while shifting it at a constant cycle.
  • the signal section cut out in this way is called a frame
  • the length of the cut-out waveform is called a frame length
  • the period of shifting the frame is called a frame period.
  • the cross-correlation function calculation unit 201 calculates a cross-correlation function of the extracted signal. Such clipping is performed on all input signals.
  • (tau) represents a lag sample and T represents a frame length.
  • the cross-correlation function calculation unit 201 calculates a cross-correlation function for all combinations of input signals. Therefore, when there are M input terminals, the number of types of cross-correlation functions is M (M ⁇ 1) / 2. For example, when there are eight input terminals, 28 types of cross-correlation functions are calculated. For the purpose of reducing the amount of calculation, it is possible to obtain a cross-correlation function for only some combinations, but the accuracy of sound source direction estimation is reduced.
  • the envelope function extraction unit 202 extracts the envelope (envelope function) of the cross correlation function supplied from the cross correlation function calculation unit 201. Envelope extraction is performed independently for all cross-correlation functions.
  • Envelope is obtained by removing a high frequency component from a waveform when the cross-correlation function is regarded as a waveform. That is, a function having a gentle shape from which a rapid change in function value is removed is defined as “correlation function envelope”. Therefore, an envelope cannot be obtained by constant multiplication or simple linear transformation.
  • the envelope extraction method there is a method of moving average the absolute value of the cross-correlation function in the lag sample axis direction.
  • the difference between the peak and the non-peak of the correlation function is significantly reduced. Therefore, the effectiveness of the method of the present embodiment in which the direction is estimated based on the peak position of the correlation function is low. Therefore, it is desirable for the envelope extraction that the difference between the peak and the non-peak is not small, but the function shape changes gradually.
  • the function zm1m2 ( ⁇ , n) obtained by taking the absolute value of the complex number of the inverse Fourier transform of the mapping function Zm1m2 (k, n) is the envelope function of the cross-correlation function ym1m2 ( ⁇ , n).
  • the envelope can also be extracted by a method using cepstrum transformation or linear predictive analysis described in Non-Patent Document 2. These methods are used in the field of audio signal processing when extracting an envelope component of a spectrum from the power spectrum of the audio signal. Compared to the Hilbert transform, the amount of calculation increases, but the envelope can be extracted with high accuracy.
  • the relative delay time calculation unit 206 obtains the relative delay time between the two microphone pairs from the input microphone position information and the sound source search target direction.
  • the relative delay time is a difference in arrival time of sound waves that is uniquely determined based on the microphone interval and the sound source direction. Assuming that the speed of sound is c, the distance between two microphones is d, and the direction of the sound source, that is, the direction of sound arrival is ⁇ , the relative delay time ⁇ , which is the difference in arrival time of sound waves, is calculated by the following equation.
  • the relative delay time calculation unit 206 calculates ⁇ for all directions and microphone pairs.
  • the direction parameter is one-dimensional.
  • the direction parameter is two-dimensional.
  • the horizontal angle is referred to as an azimuth
  • the vertical angle is referred to as an elevation angle.
  • the relative delay time calculation unit 206 calculates relative delay times corresponding to all directions for each microphone pair.
  • the search range in the horizontal direction is 0 to 90 degrees in increments of 10 degrees, that is, 0 degrees, 10 degrees, 20 degrees,... 90 degrees
  • the search range in the vertical direction is from 0 degrees in increments of 10 degrees.
  • the relative delay time calculation unit 206 calculates 70 types of relative delay times for a specific microphone pair. Since this is calculated for all microphone pairs, when the number of microphones is M, a total of 70 ⁇ M (M ⁇ 1) / 2 types of relative delay times are calculated.
  • the microphone position information input to the relative delay time calculation unit 206 is necessary to obtain the microphone interval used for calculating the relative delay time. Therefore, the microphone position information may be microphone intervals (relative positional relationship) for all microphone pairs or may be microphone position coordinates (absolute positional relationship). In the present embodiment, it is assumed that the microphone is arranged in a three-dimensional space, and the case where the microphone position information is three-dimensional coordinates will be described.
  • microphone position coordinates p m T in a three-dimensional space and a signal propagation vector u ( ⁇ , ⁇ ), which is a unit vector representing the sound source direction, are respectively defined as follows.
  • (theta) represents an azimuth and (PHI) represents an elevation angle.
  • a T represents a transposed matrix of the matrix A.
  • the integrated envelope function calculation unit 203 uses the relative delay time supplied from the relative delay time calculation unit 206 to integrate the envelope function supplied from the envelope function extraction unit 202 to calculate one integrated envelope function.
  • the integrated envelope function is generated by mixing or overlapping a plurality of envelope functions.
  • the envelope functions h m1m2 ( ⁇ , n) for all microphone pairs are superimposed.
  • the integrated envelope function for the azimuth angle ⁇ and the elevation angle ⁇ is obtained.
  • the envelope function supplied from the envelope function extraction unit 202 is z m1m2 ( ⁇ , n)
  • the integrated envelope function H n ( ⁇ , n) is calculated by the following equation.
  • the summation of sums that is, the sum of the envelope functions
  • the sum of powers may be used instead of the sum.
  • the integrated envelope function is calculated as follows.
  • the estimated direction information generation unit 204 outputs the calculated integrated envelope function as estimated direction information. Since the function value is determined for each direction (combination of the azimuth angle and the elevation angle), basically, if the function value is high, it can be determined that there is a high possibility that a sound source exists in that direction.
  • Such estimated direction information is used in various forms. For example, when the function has a plurality of peaks, it is considered that there are a plurality of sound sources having each peak as an arrival direction. Therefore, not only can the direction of each sound source be estimated simultaneously, but it can also be used to estimate the number of sound sources.
  • the possibility of the sound source based on the difference between the peak and non-peak of the function. If the difference between the peak and the non-peak is large, it can be determined that the possibility of the sound source is high. At the same time, it can be determined that the reliability of the estimated direction is high. If the number of sound sources can be assumed to be one in advance, the direction with the maximum function value may be output as estimated direction information. In this case, the estimated direction information is not the integrated envelope function but the direction itself.
  • FIG. 3 is a diagram illustrating an example of a configuration of a calculation method selection table provided in the wave source direction estimation apparatus 200 according to the present embodiment.
  • the calculation method selection table 301 stores the calculation method 312 in association with the calculation object 311. For example, when the envelope function extraction method is a calculation target, there are calculation methods such as Hilbert transform, cepstrum transform, and linear prediction analysis.
  • the wave source direction estimation apparatus 200 selects a calculation method by referring to the calculation method selection table 301. Also good.
  • FIG. 4 is a block diagram showing a hardware configuration of the wave source direction estimation apparatus 200 according to the present embodiment.
  • a CPU (Central Processing Unit) 410 is a processor for arithmetic control, and implements a functional component of the wave source direction estimation apparatus 200 of FIG. 2 by executing a program.
  • a ROM (Read Only Memory) 420 stores fixed data and programs such as initial data and programs.
  • the communication control unit 430 communicates with other devices via a network. Note that the number of CPUs 410 is not limited to one, and may be a plurality of CPUs or may include a GPU (Graphics Processing Unit) for image processing.
  • the communication control unit 430 preferably includes a CPU independent of the CPU 410 and writes or reads transmission / reception data in an area of a RAM (Random Access Memory) 440.
  • the input / output interface 460 preferably has a CPU independent of the CPU 410 and writes or reads input / output data to / from the area of the RAM 440. Therefore, the CPU 410 recognizes that the data has been received or transferred to the RAM 440 and processes the data. Further, the CPU 410 prepares the processing result in the RAM 440 and leaves the subsequent transmission or transfer to the communication control unit 430, the DMAC, or the input / output interface 460.
  • DMAC Direct Memory Access Controller
  • the RAM 440 is a random access memory that the CPU 410 uses as a temporary storage work area. In the RAM 440, an area for storing data necessary for realizing the present embodiment is secured.
  • the input signal data 441 is sound signal data collected by a sound collection device such as a microphone, or signal data input and acquired by an input signal acquisition device or the like.
  • the cross-correlation function data 442 is data related to the cross-correlation function calculated by the cross-correlation function calculation unit 201.
  • the envelope function data 443 is data related to the envelope function extracted by the envelope function extraction unit 202.
  • the integrated envelope function data 444 is data related to the integrated envelope function calculated by the integrated envelope function calculation unit 203.
  • the estimated direction information data 445 is data related to the estimated direction information generated by the estimated direction information generation unit 204.
  • the calculation method 446 is, for example, data relating to an envelope function extraction method, and is data developed from the calculation method selection table 301 shown in FIG.
  • the input / output data 447 is data input / output via the input / output interface 460.
  • the transmission / reception data 448 is data transmitted / received via the communication control unit 430.
  • the RAM 440 includes an application execution area 449 for executing various application modules.
  • the storage 450 stores a database, various parameters, or the following data or programs necessary for realizing the present embodiment.
  • the storage 450 stores a calculation method selection table 301.
  • the calculation method selection table 301 is a table for managing the relationship between the calculation target and the calculation method shown in FIG.
  • the storage 450 further stores an input signal acquisition module 451, a cross-correlation function calculation module 452, an envelope function extraction module 453, and an integrated envelope function calculation module 454.
  • the storage 450 also stores an estimated direction information generation module 455 and a relative delay time calculation module 456.
  • the input signal acquisition module 451 is a module that acquires an input signal.
  • the cross-correlation function calculation module 452 is a module that calculates a cross-correlation function based on an input signal.
  • the envelope function extraction module 453 is a module that extracts an envelope function based on the cross-correlation function.
  • the integrated envelope function calculation module 454 is a module that calculates the integrated envelope function by integrating the envelope functions.
  • the estimated direction information generation module 455 is a module that generates estimated direction information of the wave source based on the integrated envelope function.
  • the relative delay time calculation module 456 is a module for calculating a relative delay time.
  • the input / output interface 460 interfaces input / output data with input / output devices.
  • a display unit 461 and an operation unit 462 are connected to the input / output interface 460.
  • a storage medium 464 may be connected to the input / output interface 460.
  • a speaker 463 that is an audio output unit, a microphone that is an audio input unit, or a GPS position determination unit may be connected.
  • the RAM 440 and the storage 450 shown in FIG. 4 do not show programs and data related to general-purpose functions and other realizable functions that the wave source direction estimation apparatus 200 has.
  • FIG. 5 is a flowchart showing a processing procedure of the wave source direction estimation apparatus 200 according to the present embodiment. This flowchart is executed by the CPU 410 in FIG. 4 using the RAM 440, and implements a functional component of the wave source direction estimation apparatus 200 in FIG.
  • step S501 the wave source direction estimation apparatus 200 acquires an input signal.
  • step S503 the wave source direction estimation apparatus 200 calculates a cross-correlation function based on the acquired input signal.
  • step S505 the wave source direction estimation apparatus 200 extracts an envelope function based on the calculated cross-correlation function.
  • step S507 the wave source direction estimation apparatus 200 calculates a relative delay time.
  • step S509 the wave source direction estimation apparatus 200 calculates an integrated envelope function using the calculated relative delay time.
  • step S511 direction estimation information is generated based on the calculated integrated envelope function.
  • the arrival direction of the target sound included in the input signal that is, the direction in which the target object exists is estimated. This is effective in estimating the direction in which the target exists in the environment where the environmental noise level is high, using the sound generated by the target as a clue.
  • environmental noise include busy streets, streets, along highways, and places where many people and cars gather.
  • the target include humans, animals, automobiles, aircraft, ships, water bikes, and drones (small drones).
  • the position of the target sound source can be specified by performing sound source direction estimation at a plurality of locations. Thereby, even in an environment with a high environmental noise level, it is possible to accurately identify the location where a gunshot or automobile collision sound is generated.
  • FIG. 6 is a block diagram for explaining the configuration of the wave source direction estimation apparatus 600 according to this embodiment.
  • the wave source direction estimation apparatus 600 according to the present embodiment is different from the second embodiment in that it includes a conversion unit 601, a noise suppression unit 602, a cross spectrum calculation unit 603, and an envelope function extraction unit 604. Since other configurations and operations are the same as those of the second embodiment, the same configurations and operations are denoted by the same reference numerals, and detailed description thereof is omitted.
  • the conversion unit 601 converts the input signal supplied from the input terminals 20 1 to 20 M and supplies the converted signal to the noise suppression unit 602.
  • the conversion is performed for the purpose of decomposing the input signal into a plurality of frequency components.
  • a case where a representative Fourier transform is used will be described.
  • M types of input signals xm (t) are input to the conversion unit 601.
  • m is an input terminal number.
  • j represents an imaginary unit and exp represents an exponential function.
  • the noise suppression unit 602 suppresses a noise component included in the input signal based on the converted signal supplied from the conversion unit 601, and supplies the noise component to the cross spectrum calculation unit 603 as a noise suppression signal.
  • the noise suppression signal Ym (k, n) in which the signal component of Xm (k, n) below the frequency B hertz is suppressed is calculated by the following equation.
  • the signal component in that frequency band may be suppressed by the same method.
  • noise suppression methods include spectral subtraction, Wiener filter, MMSE.
  • the STSA method Minimum Mean-Square-Error Short-Time Spectral Amplitude Estimator
  • the cross spectrum calculation unit 603 calculates a cross spectrum based on the noise suppression signal supplied from the noise suppression unit 602 and transmits the cross spectrum to the envelope function extraction unit 604.
  • the cross spectrum calculation unit 603 calculates the product of two noise suppression signals for all microphone pairs.
  • a cross spectrum Pm1m2 (k, n) between the noise suppression signal Ym1 (k, n) at the input terminal m1 and the noise suppression signal Ym2 (k, n) at the input terminal m2 is calculated by the following equation.
  • the cross spectrum calculation unit 603 calculates a cross spectrum for all combinations of input signals. Therefore, when there are M input terminals, the number of types of cross spectrum is M (M ⁇ 1) / 2.
  • the cross spectrum is the Fourier transform of the cross correlation function
  • the cross spectrum can be said to be a broad cross correlation function. Accordingly, if the cross spectrum is converted into a cross-correlation function by inverse Fourier transform, the envelope function extraction unit 202 described in the second embodiment can be used instead of the envelope function extraction unit 604. In this embodiment, an example in which an envelope is extracted in a cross spectrum state will be described.
  • the envelope function extraction unit 604 extracts the envelope of the cross-correlation function based on the cross spectrum supplied from the cross spectrum calculation unit 603, and transmits it to the integrated envelope function calculation unit 203 as an envelope function. Envelope extraction is performed independently for all cross spectra. Similar to the second embodiment, an example using the Hilbert transform will be described in this embodiment. First, the following map transformation is performed.
  • the function zm1m2 ( ⁇ , n) obtained by taking the absolute value of the complex number of the inverse Fourier transform of the mapping function Zm1m2 (k, n) is the envelope function.
  • the Hilbert transform is applied to the cross spectrum that is the Fourier transform of the cross-correlation function. For this reason, the Fourier transform before performing Hilbert transform like 2nd Embodiment is unnecessary. Otherwise, the envelope can be extracted by a method using cepstrum transformation or linear predictive analysis as in the second embodiment. Although these methods require a large amount of calculation compared with the Hilbert transform, it is possible to extract the envelope with high accuracy.
  • the wave source direction estimation device 600 emphasizes the target sound by suppressing the noise component included in the input signal, and then estimates the arrival direction of the target sound, that is, the direction in which the target object exists.
  • FIG. 7 is a flowchart showing a processing procedure of the wave source direction estimation apparatus 600 according to the present embodiment. This flowchart is executed by the CPU 410 of FIG. 4 using the RAM 440, and realizes a functional component of the wave source direction estimation apparatus 600 of FIG. The same steps as those in FIG. 5 are denoted by the same step numbers and description thereof is omitted.
  • step S701 the conversion unit 601 of the wave source direction estimation apparatus 600 converts the input signal and supplies the converted input signal (converted signal) to the noise suppression unit 602.
  • step S703 the noise suppression unit 602 of the wave source direction estimation apparatus 600 suppresses the noise component included in the input signal based on the supplied converted signal, and generates a noise suppression signal.
  • step S705 the cross spectrum calculation unit 603 of the wave source direction estimation apparatus 600 calculates a cross spectrum based on the noise suppression signal, and supplies the cross spectrum to the envelope function extraction unit 604.
  • step S707 the envelope function extraction unit 604 of the wave source direction estimation apparatus 600 extracts the envelope of the cross-correlation function based on the supplied cross spectrum, and supplies it to the integrated envelope function calculation unit 203 as an envelope function.
  • This embodiment is effective when estimating the direction of the target sound source in an environment where there is a noise component having a high level of power in a specific frequency band. For example, since the power of automobile driving noise is concentrated in the low frequency range, sound source direction estimation can be performed with high accuracy in environments where the driving noise level of automobiles is high, such as along highways and highways, by removing the low frequency range. Is possible.
  • FIG. 8 is a diagram for explaining the configuration of the wave source direction estimation system 800 according to the present embodiment.
  • the wave source direction estimation system 800 according to the present embodiment uses the wave source direction estimation device according to the second embodiment or the wave source direction estimation device according to the third embodiment. Therefore, the same configurations and operations as those in the second embodiment and the third embodiment are denoted by the same reference numerals, and detailed description thereof is omitted.
  • the wave source direction estimation system 800 includes microphones 80 1 to 80 M , an AD conversion unit 801, and a display unit 802.
  • the wave source direction estimation device 600 instead of the wave source direction estimation device 200.
  • a microphone Used in place of a microphone.
  • the microphones 80 1 to 80 M convert sounds around the device including sounds generated from the target object to be estimated into electric signals and transmit them to the AD conversion unit 801.
  • the medium through which the sound is transmitted is an air medium
  • the sound reaches the microphone as air vibration.
  • the microphone converts the vibration of the air that has arrived into an electrical signal.
  • the AD conversion unit 801 converts the sound electrical signal supplied from the microphones 80 1 to 80 M into a digital signal and transmits the digital signal to the input terminals 20 1 to 20 M.
  • the display unit 802 converts the estimated direction information supplied from the wave source direction estimation device 200 (600) into visualization data such as an image and displays the visualization data on a display device such as a display.
  • the most basic visualization method is a method of projecting onto a three-dimensional solid such as a sphere. At this time, the correlation value is expressed by a color density called a heat map.
  • a method of projecting onto a two-dimensional plane such as a circle instead of a solid is also effective. There is a problem that the back side is difficult to see when displayed as a sphere. If the image is displayed on the plane projected from the front or above, the blind spot is eliminated and the listability is improved.
  • the correlation value may be expressed by contour lines instead of color shading.
  • FIG. 9A is a diagram illustrating an example of an image visualized by the display unit of the wave source direction estimation system 800 according to the present embodiment, which is obtained from estimated direction information supplied from the wave source direction estimation device 200 (600). is there.
  • FIG. 9B is a diagram illustrating another example of an image visualized on the display unit of the wave source direction estimation system 800 according to the present embodiment, in which the estimated direction information obtained by integrating the cross-correlation function instead of the envelope function is an image. It has become.
  • FIG. 9A shows that the correlation value is higher as the color is black, and the correlation value is higher as the color is white in FIG. 9B.
  • the range of azimuth is 0 to 180 degrees.
  • the elevation angle ranges from 0 to 90 degrees, and the closer to the center, the greater the elevation angle. Therefore, it is an image obtained by projecting the three-dimensional data plotted on the hemisphere onto the plane from the top. Comparing FIG. 9A and FIG. 9B, the direction of the target sound source is clear in the diagram of FIG. 9A, whereas the presence of the target sound source exists in the diagram of FIG. And the direction is unclear.
  • FIG. 10 is a diagram showing a display example of the cross-correlation function and the envelope function from which the low frequency component is removed in the wave source direction estimation system according to the present embodiment.
  • the envelope extraction which is a feature of the present embodiment, is effective when the target sound has strong power in a high frequency band. The same applies to a sound in which a high frequency component is emphasized by removing a low frequency component as in the third embodiment. This will be described with reference to FIG.
  • the peak width of the correlation function becomes narrow due to the lack of the low frequency component. If the correlation function is overlapped without performing the envelope extraction and the peak width is narrow, it is difficult to emphasize the peak by the overlap because the peak position is slightly shifted due to the influence of noise or the like. In some cases, a plurality of peaks appear, and the peak of the target sound source is not emphasized.
  • ⁇ Peak position shift is affected not only by noise but also by fluctuations in sound speed. This is because the relative delay time serving as a reference for superposition is affected by the speed of sound, as shown in Expression (3).
  • the speed of sound changes according to the temperature of a medium such as air. For this reason, in an actual environment where the temperature of air from the sound source to the microphone is not uniform, the sound speed fluctuates and a relative delay time shift occurs.
  • the peak position deviation is affected by the resolution of the estimated range, that is, the angular interval between the azimuth angle and the elevation angle.
  • the method of searching in increments of 10 degrees has been described. However, if this is increased to 20 degrees and 30 degrees, the time interval of the relative delay time becomes coarse as is apparent from equation (5).
  • FIG. 11 is a diagram showing a display example of the cross-correlation function discretized at regular intervals in the wave source direction estimation system according to the present embodiment.
  • An example in which a decrease in resolution in the estimation range causes a peak position shift will be described with reference to FIG.
  • a broken line represents a cross-correlation function before discretization.
  • FIG. 12 is a diagram showing a display example of an envelope function discretized in the wave source direction estimation system according to the present embodiment. As shown in FIG. 12, if the envelope function is discretized, a significant shift in peak position can be avoided even at the same discretization interval as in FIG.
  • the user can visually grasp the direction estimation information of the wave source.
  • the present invention may be applied to a system composed of a plurality of devices, or may be applied to a single device. Furthermore, the present invention can also be applied to a case where an information processing program that implements the functions of the embodiments is supplied directly or remotely to a system or apparatus. Therefore, in order to realize the functions of the present invention on a computer, a program installed on the computer, a medium storing the program, and a WWW (World Wide Web) server that downloads the program are also included in the scope of the present invention. . In particular, at least a non-transitory computer readable medium storing a program for causing a computer to execute the processing steps included in the above-described embodiments is included in the scope of the present invention.
  • a part or all of the above-described embodiment can be described as in the following supplementary notes, but is not limited thereto.
  • (Appendix 1) A plurality of input signal acquisition means for acquiring signals generated at the wave source as input signals; Cross-correlation function calculating means for calculating a cross-correlation function based on the input signal acquired by the input signal acquiring means; An envelope function extracting means for extracting an envelope function based on the calculated cross-correlation function; An integrated envelope function calculating means for calculating an integrated envelope function by integrating the extracted envelope functions; Based on the calculated integrated envelope function, estimated direction information generating means for generating estimated direction information of the wave source;
  • a wave source direction estimation apparatus comprising: (Appendix 2)
  • the cross-correlation function calculation means includes conversion means for converting the input signal to obtain a frequency domain signal, The wave source direction estimation apparatus according to attachment 1, wherein a cross-correlation function is calculated based on the frequency domain signal.
  • the cross-correlation function calculation means further comprises noise suppression means for obtaining an emphasized signal by suppressing noise included in the input signal based on the frequency domain signal, The wave source direction estimation apparatus according to attachment 2, wherein a cross-correlation function is calculated based on the enhancement signal.
  • the wave source direction estimation apparatus according to appendix 3 wherein the noise suppression unit suppresses noise in a predetermined frequency band.
  • a relative delay time calculating means for calculating a relative delay time of the input signal based on position information of the input signal acquiring means; The wave source direction estimation apparatus according to any one of appendices 1 to 4, wherein the integrated envelope function calculating unit calculates an integrated envelope function using the relative delay time.

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

周囲のノイズレベルが高い環境の場合であっても、波源の方向を高精度に推定することが。波源方向推定装置であって、波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、を備えた。

Description

波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム
　本発明は、波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラムに関する。
　上記技術分野において、非特許文献１には、２つ以上のマイクロフォンの受音信号を用いて複数の相互相関関数を求めてから音源（音波の発生源、発生場所）の方向を推定する技術が記載されている。具体的には、全てのマイクロフォンの組み合わせに対して、受音信号間の相互相関関数をそれぞれ求め、各組み合わせに対応する相互相関関数を総和して一つの相互相関関数を計算する。そして、総和により得られた相互相関関数が最大値を与える時間差を音波の到来時間差として算出することで、音波の到来方向を推定する。
２００１年、「マイクロホン・アレイズ」、第8章、シュプリンガー、ベルリン　ハイデルベルグ　ニューヨーク(CH.5, MICROPHONE ARRAYS, SPRINGER, BERLIN HEIDELBERG NEW YORK, 2001.) 古井、「新音響・音声工学」、近代科学社、２００６年
　しかしながら、上記文献に記載の技術では、周囲のノイズレベルが高い環境の場合には、波源の方向を高精度に推定することができなかった。
　本発明の目的は、上述の課題を解決する技術を提供することにある。
　上記目的を達成するため、本発明に係る波源方向推定装置は、
　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
　前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
　を備えた。
　上記目的を達成するため、本発明に係る波源方向推定システムは、
　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
　複数の前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
　前記推定方向情報を表示する表示手段と、
　を備えた。
　上記目的を達成するため、本発明に係る波源方向推定方法は、
　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
　を含む。
　上記目的を達成するため、本発明に係る波源方向推定プログラムは、
　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
　をコンピュータに実行させる。
　本発明によれば、周囲のノイズレベルが高い環境の場合であっても、波源の方向を高精度に推定することができる。
本発明の第１実施形態に係る波源方向推定装置の構成を示すブロック図である。 本発明の第２実施形態に係る波源方向推定装置の構成を示すブロック図である。 本発明の第２実施形態に係る波源方向推定装置の備える計算法選択テーブルの構成の一例を示す図である。 本発明の第２実施形態に係る波源方向推定装置のハードウェア構成を示すブロック図である。 本発明の第２実施形態に係る波源方向推定装置の処理手順を説明するフローチャートである。 本発明の第３実施形態に係る波源方向推定装置の構成を示すブロック図である。 本発明の第３実施形態に係る波源方向推定装置の処理手順を説明するフローチャートである。 本発明の第４実施形態に係る波源方向推定システムの構成を示すブロック図である。 本発明の第４実施形態に係る波源方向推定システムの表示部で可視化した画像の一例を示す図である。 本発明の第４実施形態に係る波源方向推定システムの表示部で可視化した画像の他の例を示す図である。 本発明の第４実施形態に係る波源方向推定システムにおいて低域成分を除去した相互相関関数および包絡関数の表示例を示す図である。 本発明の第４実施形態に係る波源方向推定システムにおいて一定間隔で離散化された相互相関関数の表示例を示す図である。 本発明の第４実施形態に係る波源方向推定システムにおいて離散化された包絡関数の表示例を示す図である。
　以下に、図面を参照して、本発明の実施の形態について例示的に詳しく説明する。ただし、以下の実施の形態に記載されている構成要素はあくまで例示であり、本発明の技術範囲をそれらのみに限定する趣旨のものではない。なお、以下の説明中における「音声信号」とは、音声その他の音響に従って生ずる直接的の電気的変化であって、音声その他の音響を伝送するためのものをいい、音声に限定されない。
　また、以下の実施形態に係る波源方向推定装置の推定対象は、空気や水の振動波である音波の発生源に限定されない。地震や地滑りなどの土や固体を媒質とする振動波の発生源にも適用できる。その場合、振動波を電気信号に変換する装置には、マイクロフォンではなく振動センサが用いられる。さらに、気体・液体・固体の振動波だけでなく、電波を用いて方向を推定する場合にも以下の実施形態に係る波源方向推定装置を適用できる。その場合、電波を電気信号に変換する装置にはアンテナが用いられる。以下の実施形態においては、波源は音源と仮定して説明する。
　［第１実施形態］
　本発明の第１実施形態としての波源方向推定装置１００について、図１を用いて説明する。波源方向推定装置１００は、取得した信号に基づいて、波源位置の方向を推定する装置である。
　図１に示すように、波源方向推定装置１００は、入力信号取得部１０１と、相互相関関数計算部１０２と、包絡関数抽出部１０３と、統合包絡関数計算部１０４と、推定方向情報生成部１０５とを含む。
　複数の入力信号取得部１０１は、波源で発生した信号を入力信号として取得する。相互相関関数計算部１０２は、複数の入力信号取得部１０１で取得した入力信号に基づいて、相互相関関数を計算する。包絡関数抽出部１０３は、計算した相互相関関数に基づいて、包絡関数を抽出する。統合包絡関数計算部１０４は、抽出した包絡関数を統合して統合包絡関数を計算する。推定方向情報生成部１０５は、計算した統合包絡関数に基づいて、波源の推定方向情報を生成する。
　本実施形態によれば、周囲のノイズレベルが高い環境の場合であっても、波源の方向を高精度に推定することができる。
　［第２実施形態］
　次に本発明の第２実施形態に係る波源方向推定装置について、図２乃至図５を用いて説明する。
　＜前提技術＞
　上記非特許文献１および非特許文献２に記載の技術では、屋外などの周囲騒音レベルが高い環境において、遠方に存在する音源の方向を高精度に推定することが困難であった。例えば、推定対象の音源（目標音源）が、マイクロフォンから遠く離れた場所に存在する場合、目標音源から放射される音の音量は、マイクロフォンに到達する時点で大幅に小さくなる。このため、目標音源の音が周囲環境雑音に埋もれて、目標音源の方向推定精度が低下することがあった。
　＜本実施形態の技術＞
　図２は、本実施形態に係る波源方向推定装置２００の構成を説明するためのブロック図である。本実施形態の波源方向推定装置２００は、例えば、デジタルビデオカメラ、スマートフォン、携帯電話、ノートパソコン、パッシブソーナーなどといった装置の一部として機能する。また、不審ドローン検知、悲鳴検知、車両事故検知などの声や音に基づいて異常を検知する異常音検知装置にも搭載される。しかし、本実施形態に係る波源方向推定装置２００の適用例はこれらに限定されるものではなく、受信音から目標音源の方向推定を要求されるあらゆる波源方向推定装置に適用可能である。
　波源方向推定装置２００は、相互相関関数計算部２０１と、包絡関数抽出部２０２と、統合包絡関数計算部２０３と、推定方向情報生成部２０４と、相対遅延時間計算部２０６とを備える。
　入力端子２０1～２０Mには、目標音源の音と集音装置であるマイクロフォン(以下、マイク)の周囲で生じる様々な雑音が混在した音とがデジタル信号(サンプル値系列)として入力される。入力端子２０1～２０Mに入力された音信号を本実施形態では入力信号と呼ぶ。そして、時刻ｔにおける入力端子２０a(ただし、ａ＝ｍ１，・・・，ｍM)の入力信号をｘａ（ｔ）と表す。
　入力端子に入力される音は、集音装置であるマイクで集音される。入力端子は複数存在するので、目標音源の音を集音する場合には、端子数と同じＭ個のマイクが同時に使用される。本実施形態では、入力端子とマイクとは一対一に対応しており、ａ番目のマイクが集音した音は、ａ番目の入力端子に供給されると仮定する。したがって、ａ番目の入力端子に供給された入力信号のことを「ａ番目のマイクの入力信号」とも呼ぶ。
　波源方向推定装置２００への入力音を集音するＭ個のマイクは、様々なレイアウトで配置される。最も基本的な配置は、一直線上に一列に並ぶ配置である。また、平面上に円状や格子状に配置する場合や、不規則な間隔でランダムに配置することもある。その他には、正多面体の頂点や球体上に配置することもある。入力端子に供給される入力信号は、このように配置されたマイクで集音された音であると仮定する。
　波源方向推定装置２００は、目標音源の音が複数のマイクに到達する時間差を利用して音源の方向を推定する。このため、マイク間隔も重要な情報となるので、波源方向推定装置２００には、入力信号だけでなく、マイク位置情報も供給される。
　入力端子がＭ個存在する場合、相互相関関数計算部２０１には、Ｍ種類の入力信号ｘａ（ｔ）が入力される。ただし、ａは入力端子番号である。相互相関関数計算部２０１は、入力端子から供給された入力信号から、適当な長さの波形を、一定の周期でずらしながら切り出す。こうして切り出した信号区間をフレーム、切り出した波形の長さをフレーム長、フレームをずらす周期をフレーム周期と呼ぶ。そして、相互相関関数計算部２０１は、切り出された信号の相互相関関数を計算する。このような切り出しは、全ての入力信号を対象に実施される。ｎ番目のフレームの時刻をｔnとすると、入力端子ｍ１の入力信号ｘｍ１（ｔ）と、入力端子ｍ２の入力信号ｘｍ２（ｔ）との相互相関関数ｙｍ１ｍ２（τ，ｎ）は、次の式で計算される。
Figure JPOXMLDOC01-appb-M000001
 ただし、τはラグサンプル、Ｔはフレーム長を表す。
　相互相関関数計算部２０１は、入力信号の全ての組み合わせに対して相互相関関数を計算する。したがって、入力端子がＭ個存在する場合、相互相関関数の種類数はＭ（Ｍ－１）／２になる。例えば入力端子が８個の場合、２８種類の相互相関関数を計算する。計算量を削減する目的で、一部の組み合わせだけに対して相互相関関数を求めることも可能であるが、音源方向の推定精度が低下する。
　包絡関数抽出部２０２は、相互相関関数計算部２０１から供給された相互相関関数の包絡（包絡関数）を抽出する。包絡の抽出は、全ての相互相関関数に対して独立に実施される。
　包絡とは、相互相関関数を波形と捉えた場合に、その波形から高周波数成分が除去されたものである。つまり、関数値の急激な時間変化が除去された、緩やかな形状を有する関数を「相関関数の包絡」と定義する。したがって、定数倍や単純な線形変換では、包絡を得ることはできない。
　最も単純な包絡抽出方法としては、相互相関関数の絶対値をラグサンプル軸方向に移動平均する方法が挙げられる。しかし移動平均では、相関関数のピークと非ピークとの差が大幅に小さくなるので、相関関数のピーク位置に基づいて方向推定を行う本実施形態の方法では有効性は低くなる。したがって、ピークと非ピークとの差は小さくならずに、関数の形状変化が緩やかになることが、包絡の抽出では望ましい。
　以上の条件を考慮して、本実施形態では、包絡の抽出にヒルベルト変換を用いた例を説明する。まず、相互相関関数ｙｍ１ｍ２（τ，ｎ）のフーリエ変換を計算する。求めたフーリエ変換をＹｍ１ｍ２（ｋ，ｎ）（ただし、ｋは周波数ビン番号、ｋ＝０，１，・・・，Ｋ－１）とすると、以下の写像変換を行う。
Figure JPOXMLDOC01-appb-M000002
 　そして、写像関数Ｚｍ１ｍ２（ｋ，ｎ）の逆フーリエ変換の複素数の絶対値を取って得られた関数ｚｍ１ｍ２（τ，ｎ）が、相互相関関数ｙｍ１ｍ２（τ，ｎ）の包絡関数となる。
　他には、非特許文献２に記載されているケプストラム変換や、線型予測分析を用いた方法でも包絡を抽出できる。これらの方法は、音声信号処理の分野において、音声信号のパワースペクトルからスペクトルの包絡成分を抽出する際に用いられている。ヒルベルト変換と比較すると、計算量は多くなるものの、包絡を高精度に抽出することが可能である。
　相対遅延時間計算部２０６は、入力されたマイク位置情報と音源探索対象方向とから、２つのマイクペア間の相対遅延時間を求める。相対遅延時間とは、マイク間隔と音源方向とに基づいて一意に定まる音波の到達時間差のことである。音速をｃと仮定し、ある２つのマイクの間隔をｄ、音源方向、つまり音の到来方向をθとした場合、音波の到達時間差である相対遅延時間τは次の式で計算される。
Figure JPOXMLDOC01-appb-M000003
 　相対遅延時間計算部２０６では、全ての方向とマイクペアに対して、上記のτを計算する。
　音源探索方向が１軸である場合、例えば垂直方向を無視して水平方向だけを対象とする場合は、方向パラメータは一次元である。一方、２軸である場合、つまり水平方向と垂直方向との両方の角度を用いて音源方向を推定する場合、方向パラメータは二次元となる。本実施形態では、方向パラメータが２次元の場合について説明する。また、水平方向の角度を方位角、垂直方向の角度を仰角と呼ぶこととする。
　音源方向推定では、探索範囲にある全ての方向に対して、包絡関数を統合したものである統合包絡関数を計算する。そこで、相対遅延時間計算部２０６では、各マイクペアに対して、全ての方向に対応する相対遅延時間を計算する。例えば、水平方向の探索範囲が１０度刻みで０度から９０度まで、つまり０度、１０度、２０度、・・・９０度であり、垂直方向の探索範囲が１０度刻みで０度から６０度まで、つまり０度、１０度、２０度、・・・６０度である場合を考える。この場合、相対遅延時間計算部２０６は、ある特定のマイクペアに対して７０種類の相対遅延時間を計算する。これを全てのマイクペアに対して計算するので、マイク数がＭ個の場合、相対遅延時間は合計で７０・Ｍ（Ｍ－１）／２種類計算される。
　相対遅延時間計算部２０６に入力されるマイク位置情報は、相対遅延時間の計算に用いるマイク間隔を得るために必要となる。したがって、マイク位置情報は、全てのマイクペアに対するマイク間隔(相対的な位置関係)の場合もあれば、マイク位置座標(絶対的な位置関係)の場合もある。本実施形態では、マイクが３次元空間上に配置されていると仮定し、マイク位置情報が３次元座標である場合について説明する。
　相対遅延時間の計算方法は次の通りである。まず、三次元空間におけるマイク位置座標ｐｍ Ｔ、音源方向を表す単位ベクトルである信号伝搬ベクトルｕ(θ、Φ)をそれぞれ次のように定義する。
Figure JPOXMLDOC01-appb-M000004
 　ただし、θは方位角、Φは仰角を表す。また、ＡＴは、行列Ａの転置行列を表す。このとき、マイクｍ１とマイクｍ２との組み合わせに対する相対遅延時間τm1m2（θ、Φ）は、次のように与えられる。
Figure JPOXMLDOC01-appb-M000005
 　ただし、ｃは音速である。式（５）で得られる相対遅延時間は実数なので、統合包絡関数計算部２０３には整数化した相対遅延時間τm1m2（θ、Φ）が供給される。
　統合包絡関数計算部２０３は、相対遅延時間計算部２０６から供給された相対遅延時間を用いて、包絡関数抽出部２０２から供給された包絡関数を統合して、一つの統合包絡関数を算出する。統合包絡関数は、複数の包絡関数を混合したり、重ね合わせることにより生成される。相対遅延時間計算部２０６から供給された、方位角θと仰角Φとに対する相対遅延時間τm1m2（θ、Φ）に基づいて、全てのマイクペアに対する包絡関数ｈm1m2（τ、n）を重ね合わせて、方位角θと仰角Φとに対する統合包絡関数求める。包絡関数抽出部２０２から供給された包絡関数を ｚm1m2（τ、n）とすると、統合包絡関数 Ｈn（τ、n）は次式で計算される。
Figure JPOXMLDOC01-appb-M000006
 　すなわち、和による重ね合わせ、つまり各包絡関数の総和を統合包絡関数とする。また、総和では無く総乗を用いてもよい。その場合、統合包絡関数は次のように計算される。
Figure JPOXMLDOC01-appb-M000007
 　この他にも、包絡関数を個別に定数倍してから総和を求める重み付け和や、包絡関数の指数や対数を求めてから総和や総乗を計算する方法も同様に用いることが可能である。
　そして、推定方向情報生成部２０４は、算出された統合包絡関数を推定方向情報として出力する。各方向(方位角と仰角との組み合わせ)に対して関数値が定まるので、基本的に関数値が高ければ、その方向に音源が存在する可能性が高いと判断できる。
　このような推定方向情報は、様々な形で利用される。例えば、関数が複数のピークを有する場合には、各ピークを到来方向とする音源が複数存在すると考えられる。したがって、各音源の方向を同時に推定できるだけでなく、音源数の推定にも用いることが可能である。
　また、関数のピークと非ピークとの差分に基づき、音源の存在可能性を判定することも可能である。もし、ピークと非ピークとの差分が大きければ音源の存在可能性が高いと判定できる。同時に、推定方向の信頼性も高いと判断できる。もし、音源数が１つと事前に仮定できる場合には、関数値が最大となる方向を推定方向情報として出力してもよい。この場合、推定方向情報は、統合包絡関数では無く、方向そのものとなる。
　図３は、本実施形態に係る波源方向推定装置２００の備える計算法選択テーブルの構成の一例を示す図である。計算法選択テーブル３０１は、計算対象３１１に関連付けて計算法３１２を記憶する。例えば、包絡関数抽出法が計算対象の場合、ヒルベルト変換やケプストラム変換、線型予測分析などの計算法があり、波源方向推定装置２００は、計算法選択テーブル３０１を参照して計算法を選択してもよい。
　図４は、本実施形態に係る波源方向推定装置２００のハードウェア構成を示すブロック図である。
　ＣＰＵ(Central Processing Unit)４１０は演算制御用のプロセッサであり、プログラムを実行することで図２の波源方向推定装置２００の機能構成部を実現する。ＲＯＭ(Read Only Memory)４２０は、初期データおよびプログラムなどの固定データおよびプログラムを記憶する。また、通信制御部４３０は、ネットワークを介してその他の装置などと通信する。なお、ＣＰＵ４１０は１つに限定されず、複数のＣＰＵであっても、あるいは画像処理用のＧＰＵ(Graphics Processing Unit)を含んでもよい。また、通信制御部４３０は、ＣＰＵ４１０とは独立したＣＰＵを有して、ＲＡＭ(Random Access Memory)４４０の領域に送受信データを書き込みあるいは読み出しするのが望ましい。また、ＲＡＭ４４０とストレージ４５０との間でデータを転送するＤＭＡＣ(Direct Memory Access Controller)を設けるのが望ましい（図示なし）。さらに、入出力インタフェース４６０は、ＣＰＵ４１０とは独立したＣＰＵを有して、ＲＡＭ４４０の領域に入出力データを書き込みあるいは読み出しするのが望ましい。したがって、ＣＰＵ４１０は、ＲＡＭ４４０にデータが受信あるいは転送されたことを認識してデータを処理する。また、ＣＰＵ４１０は、処理結果をＲＡＭ４４０に準備し、後の送信あるいは転送は通信制御部４３０やＤＭＡＣ、あるいは入出力インタフェース４６０に任せる。
　ＲＡＭ４４０は、ＣＰＵ４１０が一時記憶のワークエリアとして使用するランダムアクセスメモリである。ＲＡＭ４４０には、本実施形態の実現に必要なデータを記憶する領域が確保されている。入力信号データ４４１は、マイクなどの集音装置で集音した音信号データや入力信号取得装置などに入力され、取得された信号データなどである。
　相互相関関数データ４４２は、相互相関関数計算部２０１で計算した相互相関関数に関するデータである。包絡関数データ４４３は、包絡関数抽出部２０２で抽出した包絡関数に関するデータである。統合包絡関数データ４４４は、統合包絡関数計算部２０３で計算した統合包絡関数に関するデータである。推定方向情報データ４４５は、推定方向情報生成部２０４で生成した推定方向情報に関するデータである。計算法４４６は、例えば、包絡関数の抽出法に関するデータであり、図３に示した計算法選択テーブル３０１から展開されたデータである。
　入出力データ４４７は、入出力インタフェース４６０を介して入出力されるデータである。送受信データ４４８は、通信制御部４３０を介して送受信されるデータである。また、ＲＡＭ４４０は、各種アプリケーションモジュールを実行するためのアプリケーション実行領域４４９を有する。
　ストレージ４５０には、データベースや各種のパラメータ、あるいは本実施形態の実現に必要な以下のデータまたはプログラムが記憶されている。ストレージ４５０は、計算法選択テーブル３０１を格納する。計算法選択テーブル３０１は、図３に示した、計算対象と計算法との関係を管理するテーブルである。
　ストレージ４５０は、さらに、入力信号取得モジュール４５１と、相互相関関数計算モジュール４５２と、包絡関数抽出モジュール４５３と、統合包絡関数計算モジュール４５４とを格納する。また、ストレージ４５０は、推定方向情報生成モジュール４５５と、相対遅延時間計算モジュール４５６とを格納する。
　入力信号取得モジュール４５１は、入力信号を取得するモジュールである。相互相関関数計算モジュール４５２は、入力信号に基づいて相互相関関数を計算するモジュールである。包絡関数抽出モジュール４５３は、相互相関関数に基づいて包絡関数を抽出するモジュールである。統合包絡関数計算モジュール４５４は、包絡関数を統合して統合包絡関数を計算するモジュールである。推定方向情報生成モジュール４５５は、統合包絡関数に基づいて波源の推定方向情報を生成するモジュールである。相対遅延時間計算モジュール４５６は、相対遅延時間を計算するモジュールである。これらのモジュール４５１～４５６は、ＣＰＵ４１０によりＲＡＭ４４０のアプリケーション実行領域４４９に読み出され、実行される。制御プログラム４５７は、波源方向推定装置２００の全体を制御するためのプログラムである。
　入出力インタフェース４６０は、入出力機器との入出力データをインタフェースする。入出力インタフェース４６０には、表示部４６１、操作部４６２、が接続される。また、入出力インタフェース４６０には、さらに、記憶媒体４６４が接続されてもよい。さらに、音声出力部であるスピーカ４６３や、音声入力部であるマイク、あるいは、ＧＰＳ位置判定部が接続されてもよい。なお、図４に示したＲＡＭ４４０やストレージ４５０には、波源方向推定装置２００が有する汎用の機能や他の実現可能な機能に関するプログラムやデータは図示されていない。
　図５は、本実施形態に係る波源方向推定装置２００の処理手順を示すフローチャートである。このフローチャートは、図４のＣＰＵ４１０がＲＡＭ４４０を使用して実行し、図２の波源方向推定装置２００の機能構成部を実現する。
　ステップＳ５０１において、波源方向推定装置２００は、入力信号を取得する。ステップＳ５０３において、波源方向推定装置２００は、取得した入力信号に基づいて相互相関関数を計算する。ステップＳ５０５において、波源方向推定装置２００は、計算した相互相関関数に基づいて、包絡関数を抽出する。ステップＳ５０７において、波源方向推定装置２００は、相対遅延時間を計算する。ステップＳ５０９において、波源方向推定装置２００は、計算した相対遅延時間を用いて統合包絡関数を計算する。ステップＳ５１１において、計算した統合包絡関数に基づいて、方向推定情報を生成する。
　本実施形態によれば、入力信号に含まれる目標音の到来方向、すなわち目標物体が存在する方向を推定する。環境雑音レベルが高い環境において、目標物が発する音を手掛かりに、目標物が存在する方向を推定する場合に有効である。環境雑音の例としては、繁華街や街頭、街道沿い、人や自動車が多く集まる場所が挙げられる。また、目標物の例としては、人間や動物、自動車、航空機、船舶、水上バイク、ドローン（小型無人機）が挙げられる。
　例えば、屋外のテーマパークや展示会場などに接近する不審な自動車・船舶・ドローンなどを検知し、それらの方向を推定することで、不審者や不審物を効率的に取り締ることが可能である。また、音源方向推定を複数箇所で実施することで、目標音源の位置を特定できる。これにより、環境雑音レベルが高い環境でも、銃声や自動車の衝突音の発生箇所などを正確に特定することが可能となる。
　［第３実施形態］
　次に本発明の第３実施形態に係る波源方向推定装置について、図６および図７を用いて説明する。図６は、本実施形態に係る波源方向推定装置６００の構成を説明するためのブロック図である。本実施形態に係る波源方向推定装置６００は、上記第２実施形態と比べると、変換部６０１、雑音抑圧部６０２、クロススペクトル計算部６０３および包絡関数抽出部６０４を有する点で異なる。その他の構成および動作は、第２実施形態と同様であるため、同じ構成および動作については同じ符号を付してその詳しい説明を省略する。
　変換部６０１は、入力端子２０1～２０Mから供給された入力信号を変換し、雑音抑圧部６０２へ供給する。変換は、入力信号を複数の周波数成分に分解する目的で実施される。ここでは、代表的なフーリエ変換を用いた場合について説明する。
　第２実施形態の相互相関関数計算部２０１と同様に、変換部６０１にはＭ種類の入力信号ｘｍ（ｔ）が入力される。ただし、ｍは入力端子番号である。変換部６０１は、入力端子から供給された入力信号から、適当な長さの波形を一定の周期でずらしながら切り出す。そして、フーリエ変換を用いて切り出された信号を周波数領域信号に変換する。ｎをフレーム番号とし、切り出した入力信号ｘｍ（ｔ，ｎ）（ｔ＝０，１，・・・，Ｔ－１）とすると、ｘｍ（ｔ，ｎ）のフーリエ変換Ｘｍ（ｋ，ｎ）は次のように計算される。
Figure JPOXMLDOC01-appb-M000008
 　ただし、ｊは虚数単位、ｅｘｐは指数関数を表す。
　雑音抑圧部６０２は、変換部６０１から供給された変換信号に基づいて、入力信号に含まれる雑音成分を抑圧し、雑音抑圧信号としてクロススペクトル計算部６０３へ供給する。
　雑音成分が低周波数帯域に存在していることが事前に判明している場合には、低周波数の信号成分をゼロにする方法が有効である。入力信号のサンプリングレートをＦｓとすると、Ｘｍ（ｋ，ｎ）の周波数Ｂヘルツ以下の信号成分が抑圧された雑音抑圧信号Ｙｍ（ｋ，ｎ）は次の式で計算される。
Figure JPOXMLDOC01-appb-M000009
 　もし、特定の周波数帯域に雑音が存在する場合には、その周波数帯域の信号成分を同様の方法で抑圧すればよい。
　この他には、入力信号に目的音源の音が含まれていない時間帯に、入力信号の平均パワースペクトルを求め、その平均パワースペクトルに基づいて雑音成分を抑圧する方法がある。このような雑音抑圧方法には、スペクトル減算法、ウィーナーフィルタ法、ＭＭＳＥ
　ＳＴＳＡ法（Minimum Mean-Square-Error Short-Time Spectral Amplitude Estimator）などが知られている。マイク設置場所に定常的に存在する環境雑音を抑圧する場合に有効である。
　クロススペクトル計算部６０３は、雑音抑圧部６０２から供給された雑音抑圧信号に基づいて、クロススペクトルを計算し、包絡関数抽出部６０４へ伝達する。クロススペクトル計算部６０３は、全てのマイクペアに対して２つの雑音抑圧信号の積を計算する。入力端子ｍ１の雑音抑圧信号Ｙｍ１（ｋ，ｎ）と、入力端子ｍ２の雑音抑圧信号Ｙｍ２（ｋ，ｎ）とのクロススペクトルＰｍ１ｍ２（ｋ，ｎ）は、次の式で計算される。
Figure JPOXMLDOC01-appb-M000010
 　ただし、Ｘバーは、Ｘの複素共役を表す。クロススペクトル計算部６０３は、入力信号の全ての組み合わせに対してクロススペクトルを計算する。したがって、入力端子がＭ個存在する場合、クロススペクトルの種類数はＭ（Ｍ－１）／２になる。
　クロススペクトルは相互相関関数のフーリエ変換であるので、クロススペクトルも広義の相互相関関数と言える。したがって、ここでクロススペクトルを逆フーリエ変換により相互相関関数に変換すれば、包絡関数抽出部６０４の代わりに、第２実施形態に記載した包絡関数抽出部２０２を用いることが可能である。本実施形態では、クロススペクトルの状態で包絡を抽出する例について説明する。
　包絡関数抽出部６０４は、クロススペクトル計算部６０３から供給されたクロススペクトルに基づいて、相互相関関数の包絡を抽出し、包絡関数として統合包絡関数計算部２０３へ伝達する。包絡の抽出は、全てのクロススペクトルに対して独立に実施される。第２実施形態と同様に、本実施形態でもヒルベルト変換を用いた例を説明する。はじめに、以下の写像変換を行う。
Figure JPOXMLDOC01-appb-M000011
 　そして、写像関数Ｚｍ１ｍ２（ｋ，ｎ）の逆フーリエ変換の複素数の絶対値を取って得られた関数ｚｍ１ｍ２（τ，ｎ）が、包絡関数となる。本実施形態では、相互相関関数のフーリエ変換であるクロススペクトルに対してヒルベルト変換を適用している。このため、第２実施形態のようにヒルベルト変換を行う前のフーリエ変換が不要である。他には、第２実施形態と同様にケプストラム変換や、線型予測分析を用いた方法でも包絡を抽出できる。これらの方法は、ヒルベルト変換と比較して計算量が多くなるものの、包絡を高精度に抽出することが可能である。
　以上の構成により、波源方向推定装置６００は、入力信号に含まれる雑音成分の抑圧により目標音を強調してから、目標音の到来方向、すなわち目標物体が存在する方向を推定する。
　図７は、本実施形態に係る波源方向推定装置６００の処理手順を示すフローチャートである。このフローチャートは、図４のＣＰＵ４１０がＲＡＭ４４０を使用して実行し、図６の波源方向推定装置６００の機能構成部を実現する。なお、図５と同様の処理には同じステップ番号を付して説明を省略する。
　ステップＳ７０１において、波源方向推定装置６００の変換部６０１は、入力信号を変換して、変換した入力信号（変換信号）を雑音抑圧部６０２へ供給する。ステップＳ７０３において、波源方向推定装置６００の雑音抑圧部６０２は、供給された変換信号に基づいて、入力信号に含まれる雑音成分を抑圧して、雑音抑圧信号を生成する。ステップＳ７０５において、波源方向推定装置６００のクロススペクトル計算部６０３は、雑音抑圧信号に基づいて、クロススペクトルを計算し、包絡関数抽出部６０４へ供給する。ステップＳ７０７において、波源方向推定装置６００の包絡関数抽出部６０４は、供給されたクロススペクトルに基づいて、相互相関関数の包絡を抽出し、包絡関数として統合包絡関数計算部２０３へ供給する。
　本実施形態によれば、特定の周波数帯域に高いレベルのパワーを有する雑音成分が存在する環境において、目標音源の方向を推定する場合に有効である。例えば、自動車走行雑音のパワーは低域に集中していることから、低域を除去することで、高速道路や街道沿いなどの自動車の走行雑音レベルが高い環境で、高精度に音源方向推定を行うことが可能である。
　また、目標音のレベルが環境雑音よりも高い帯域が事前に判明している場合にも有効である。例えば、バイクやヘリコプター、ドローンは人間の声や自動車走行音と比べると、高周波数帯域のパワーが大きい。このため、バイクやヘリコプターなどの方向を推定したい場合にも、低域成分を除去し、目標音が支配的な高域成分を強調することで、音源方向推定の推定精度が向上する。
　［第４実施形態］
　次に本発明の第４実施形態に係る波源方向推定システムについて、図８乃至図１２を用いて説明する。図８は、本実施形態に係る波源方向推定システム８００の構成を説明するための図である。本実施形態に係る波源方向推定システム８００は、上記第２実施形態に係る波源方向推定装置または第３実施形態に係る波源方向推定装置を用いている。したがって、第２実施形態および第３実施形態と同じ構成および動作については同じ符号を付してその詳しい説明を省略する。
　本実施形態に係る波源方向推定システム８００は、マイク８０1～８０M、ＡＤ変換部８０１および表示部８０２を有する。なお、本実施形態では、波源方向推定装置２００の代わりに波源方向推定装置６００を用いることが可能である。また、波源が音源であるという仮定で説明するため、マイクを用いた例を説明するが、音源以外の場合には、その波源から放射される波動を受信して電気信号に変換できる各種センサが、マイクの代わりに用いられる。
　マイク８０1～８０Mは、推定対象である目標物体から生じる音を含めた装置周辺の音を電気信号に変換し、ＡＤ変換部８０１へ伝達する。音が伝わる媒質が空気媒質である場合、音は空気の振動としてマイクに到達する。マイクは、到達した空気の振動を電気信号に変換する。
　ＡＤ変換部８０１は、マイク８０1～８０Mから供給された音の電気信号をデジタル信号に変換し、入力端子２０1～２０Mに伝達する。
　表示部８０２は、波源方向推定装置２００（６００）から供給された推定方向情報を画像などの可視化データに変換し、ディスプレイなどの表示装置に表示する。最も基本的な可視化方法は、球体などの３次元の立体に投影する方法である。その際、相関値をヒートマップと呼ばれる色の濃淡などで表現する。立体では無く、円などの２次元平面に投影する方法も有効である。球体だと表示されたときに裏側が見づらいという問題がある。正面や上から投影した平面上に表示すれば、死角が無くなり一覧性が向上する。相関値は、色の濃淡ではなく、等高線で表現してもよい。相関値の時間変化を把握しやすくする上では、横軸を方向、縦軸を時間にした二次元平面にプロットした方法も有効である。その際、仰角と方位角は別々にプロットしたほうがよい。
　図９Ａは、本実施形態に係る波源方向推定システム８００の表示部で可視化した画像の一例を示す図であり、波源方向推定装置２００（６００）から供給された推定方向情報から得られたものである。図９Ｂは、本実施形態に係る波源方向推定システム８００の表示部で可視化した画像の他の例を示す図であり、包絡関数ではなく、相互相関関数を統合して得た推定方向情報を画像化したものである。
　これらは、本実施形態の特徴である包絡抽出の効果を確認する目的で、同一の入力信号から同一のフレーム時刻に取得した。例の作成には、正面方向から接近するドローン(小型無人機)の飛翔音を用いた。多面体の各頂点に配置された複数のマイクを用いて飛翔音を集音した。
　図９Ａは、色が黒いほど相関値が高いことを表し、図９Ｂでは色が白いほど相関値が高くなっている。方位角の範囲は０から１８０度である。また、仰角の範囲は０から９０度で、中心に近いほど仰角が大きくなる。したがって、半球体上にプロットした３次元データを上から平面上に投影した画像となっている。図９Ａと図９Ｂとを比較すると、図９Ａの図は、目標音源の方向が明確であるのに対して、図９Ｂの図は、同等の白さの箇所が複数あるため、目標音源の存在とその方向が不明確であることが分かる。
　図１０は、本実施形態に係る波源方向推定システムにおいて低域成分を除去した相互相関関数および包絡関数の表示例を示す図である。本実施形態の特徴である包絡の抽出は、目標音が高い周波数帯域に強いパワーを有する場合に有効である。このことは、第３実施形態のように低域成分の除去により、高域成分を強調した音の場合でも同様である。これを図１０を用いて説明する。
　図１０に示す通り、高周波数成分のみが存在する音の相互相関関数は、低域成分の欠落が原因で、相関関数のピーク幅が狭くなる。もし包絡抽出を行わず、ピーク幅が狭い状態で相関関数を重ね合わせると、雑音などの影響によりピーク位置が僅かにずれるだけで、重ね合わせによるピークの強調が困難になる。場合によっては、複数のピークが出現し、目標音源のピークが強調されなくなる。
　ピーク位置のずれは、雑音だけでなく音速のゆらぎにも影響を受ける。式（３）に示す通り、重ね合わせの基準となる相対遅延時間は、音速の影響を受けるためである。音速は、空気などの媒質の温度に応じて変化する。このため、音源からマイクまでの空気の温度が一様ではない実環境では、音速に揺らぎが生じ、相対遅延時間のずれがもたらされる。
　また、ピーク位置のずれは、推定範囲の解像度、つまり方位角と仰角との角度間隔にも影響を受ける。第２実施形態では１０度刻みに探索する方法を説明したが、これを２０度、３０度と大きくすると、式（５）から明らかなように相対遅延時間の時間間隔が粗くなる。
　図１１は、本実施形態に係る波源方向推定システムにおいて一定間隔で離散化された相互相関関数の表示例を示す図である。推定範囲の解像度低下がピーク位置のずれをもたらす例について、図１１を用いて説明する。破線は、離散化前の相互相関関数を表す。このように離散化の間隔が大きいと、相関関数のピーク位置は、本来とは大きく異なる位置に一致することがある。このずれは離散化の間隔が小さくなるにつれて発生しにくくなるが、相関関数の重ね合わせに要する計算量が大幅に増大する。
　図１２は、本実施形態に係る波源方向推定システムにおいて離散化された包絡関数の表示例を示す図である。図１２に示すように、包絡関数を離散化するのであれば、図１１と同じ離散化の間隔であってもピーク位置の大幅なずれを回避できる。
　本実施形態によれば、推定方向情報を画像などの可視化データとして表示するので、ユーザが波源の方向推定情報を視覚的に把握することができる。
　［他の実施形態］
　以上、実施形態を参照して本願発明を説明したが、本願発明は上記実施形態に限定されるものではない。本願発明の構成や詳細には、本願発明のスコープ内で当業者が理解し得る様々な変更をすることができる。また、それぞれの実施形態に含まれる別々の特徴を如何様に組み合わせたシステムまたは装置も、本発明の範疇に含まれる。
　また、本発明は、複数の機器から構成されるシステムに適用されてもよいし、単体の装置に適用されてもよい。さらに、本発明は、実施形態の機能を実現する情報処理プログラムが、システムあるいは装置に直接あるいは遠隔から供給される場合にも適用可能である。したがって、本発明の機能をコンピュータで実現するために、コンピュータにインストールされるプログラム、あるいはそのプログラムを格納した媒体、そのプログラムをダウンロードさせるＷＷＷ(World Wide Web)サーバも、本発明の範疇に含まれる。特に、少なくとも、上述した実施形態に含まれる処理ステップをコンピュータに実行させるプログラムを格納した非一時的コンピュータ可読媒体（non-transitory computer readable medium）は本発明の範疇に含まれる。
　［実施形態の他の表現］
　上記の実施形態の一部または全部は、以下の付記のようにも記載されうるが、以下には限られない。
（付記１）
　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
　前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
　を備えた波源方向推定装置。
（付記２）
　前記相互相関関数計算手段は、前記入力信号を変換して周波数領域信号を得る変換手段を備え、
　前記周波数領域信号に基づいて、相互相関関数を計算する付記１に記載の波源方向推定装置。
（付記３）
　前記相互相関関数計算手段は、前記周波数領域信号に基づいて、前記入力信号に含まれる雑音を抑圧して強調信号を得る雑音抑圧手段をさらに備え、
　前記強調信号に基づいて、相互相関関数を計算する付記２に記載の波源方向推定装置。
（付記４）
　前記雑音抑圧手段は、所定の周波数帯域の雑音を抑圧する付記３に記載の波源方向推定装置。
（付記５）
　前記入力信号取得手段の位置情報に基づいて、前記入力信号の相対遅延時間を計算する相対遅延時間計算手段をさらに備え、
　前記統合包絡関数計算手段は、前記相対遅延時間を用いて、統合包絡関数を計算する付記１乃至４のいずれか１項に記載の波源方向推定装置。
（付記６）
　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
　複数の前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
　前記推定方向情報を表示する表示手段と、
　を備えた波源方向推定システム。
（付記７）
　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
　を含む波源方向推定方法。
（付記８）
　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
　をコンピュータに実行させる波源方向推定プログラム。
　この出願は、２０１６年６月１５日に出願された日本出願特願２０１６－１１８４５５を基礎とする優先権を主張し、その開示の全てをここに取り込む。

Claims (8)

  1. 　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
    　前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
    　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
    　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
    　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
    　を備えた波源方向推定装置。
  2. 　前記相互相関関数計算手段は、前記入力信号を変換して周波数領域信号を得る変換手段を備え、
    　前記周波数領域信号に基づいて、相互相関関数を計算する請求項１に記載の波源方向推定装置。
  3. 　前記相互相関関数計算手段は、前記周波数領域信号に基づいて、前記入力信号に含まれる雑音を抑圧して強調信号を得る雑音抑圧手段をさらに備え、
    　前記強調信号に基づいて、相互相関関数を計算する請求項２に記載の波源方向推定装置。
  4. 　前記雑音抑圧手段は、所定の周波数帯域の雑音を抑圧する請求項３に記載の波源方向推定装置。
  5. 　前記入力信号取得手段の位置情報に基づいて、前記入力信号の相対遅延時間を計算する相対遅延時間計算手段をさらに備え、
    　前記統合包絡関数計算手段は、前記相対遅延時間を用いて、統合包絡関数を計算する請求項１乃至４のいずれか１項に記載の波源方向推定装置。
  6. 　波源で発生した信号を入力信号として取得する複数の入力信号取得手段と、
    　複数の前記入力信号取得手段で取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算手段と、
    　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出手段と、
    　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算手段と、
    　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成手段と、
    　前記推定方向情報を表示する表示手段と、
    　を備えた波源方向推定システム。
  7. 　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
    　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
    　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
    　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
    　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
    　を含む波源方向推定方法。
  8. 　波源で発生した信号を入力信号として取得する入力信号取得ステップと、
    　前記入力信号取得ステップで取得した前記入力信号に基づいて、相互相関関数を計算する相互相関関数計算ステップと、
    　計算した前記相互相関関数に基づいて、包絡関数を抽出する包絡関数抽出ステップと、
    　抽出した前記包絡関数を統合して統合包絡関数を計算する統合包絡関数計算ステップと、
    　計算した前記統合包絡関数に基づいて、前記波源の推定方向情報を生成する推定方向情報生成ステップと、
    　をコンピュータに実行させる波源方向推定プログラム。
PCT/JP2017/002787 2016-06-15 2017-01-26 波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム WO2017216999A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2018523292A JP6874766B2 (ja) 2016-06-15 2017-01-26 波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム
US16/308,860 US11454694B2 (en) 2016-06-15 2017-01-26 Wave source direction estimation apparatus, wave source direction estimation system, wave source direction estimation method, and wave source direction estimation program

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2016-118455 2016-06-15
JP2016118455 2016-06-15

Publications (1)

Publication Number Publication Date
WO2017216999A1 true WO2017216999A1 (ja) 2017-12-21

Family

ID=60663497

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2017/002787 WO2017216999A1 (ja) 2016-06-15 2017-01-26 波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム

Country Status (3)

Country Link
US (1) US11454694B2 (ja)
JP (1) JP6874766B2 (ja)
WO (1) WO2017216999A1 (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111505568A (zh) * 2020-04-14 2020-08-07 哈尔滨工程大学 基于差分优化的四面体阵目标方位估计方法
JP2020176902A (ja) * 2019-04-17 2020-10-29 株式会社Ihi 時間差測定装置及び到来方向推定装置

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11408963B2 (en) * 2018-06-25 2022-08-09 Nec Corporation Wave-source-direction estimation device, wave-source-direction estimation method, and program storage medium
US11295543B2 (en) * 2020-03-31 2022-04-05 International Business Machines Corporation Object detection in an image
CN117370731B (zh) * 2023-10-10 2024-06-04 广州远动信息技术有限公司 一种基于卷积神经网络的声音到达时间估计方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH11304906A (ja) * 1998-04-20 1999-11-05 Nippon Telegr & Teleph Corp <Ntt> 音源位置推定方法およびそのプログラムを記録した記録媒体
US20030081503A1 (en) * 2001-10-23 2003-05-01 Barnard Thomas J. System and method for localizing targets using multiple arrays
JP2010185690A (ja) * 2009-02-10 2010-08-26 Mitsubishi Electric Corp 方位探知装置及び方位探知方法
WO2015137425A1 (ja) * 2014-03-14 2015-09-17 株式会社村田製作所 超音波距離測定装置および超音波距離測定方法
JP2016114512A (ja) * 2014-12-16 2016-06-23 日本電気株式会社 振動発生源推定装置、方法およびプログラム

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH11304906A (ja) * 1998-04-20 1999-11-05 Nippon Telegr & Teleph Corp <Ntt> 音源位置推定方法およびそのプログラムを記録した記録媒体
US20030081503A1 (en) * 2001-10-23 2003-05-01 Barnard Thomas J. System and method for localizing targets using multiple arrays
JP2010185690A (ja) * 2009-02-10 2010-08-26 Mitsubishi Electric Corp 方位探知装置及び方位探知方法
WO2015137425A1 (ja) * 2014-03-14 2015-09-17 株式会社村田製作所 超音波距離測定装置および超音波距離測定方法
JP2016114512A (ja) * 2014-12-16 2016-06-23 日本電気株式会社 振動発生源推定装置、方法およびプログラム

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2020176902A (ja) * 2019-04-17 2020-10-29 株式会社Ihi 時間差測定装置及び到来方向推定装置
JP7267825B2 (ja) 2019-04-17 2023-05-02 株式会社Ihi 時間差測定装置及び到来方向推定装置
CN111505568A (zh) * 2020-04-14 2020-08-07 哈尔滨工程大学 基于差分优化的四面体阵目标方位估计方法

Also Published As

Publication number Publication date
JP6874766B2 (ja) 2021-05-19
US11454694B2 (en) 2022-09-27
US20190146054A1 (en) 2019-05-16
JPWO2017216999A1 (ja) 2019-04-04

Similar Documents

Publication Publication Date Title
WO2017216999A1 (ja) 波源方向推定装置、波源方向推定システム、波源方向推定方法および波源方向推定プログラム
JP6769495B2 (ja) 相関関数生成装置、相関関数生成方法、相関関数生成プログラムおよび波源方向推定装置
JP6526083B2 (ja) 源信号分離のためのシステム及び方法
US9432770B2 (en) Method and device for localizing sound sources placed within a sound environment comprising ambient noise
US20170325023A1 (en) Multi-microphone neural network for sound recognition
CN105301563B (zh) 一种基于一致聚焦变换最小二乘法的双声源定位方法
CN104200813B (zh) 基于声源方向实时预测跟踪的动态盲信号分离方法
CN102305925A (zh) 一种机器人连续声源定位方法
CN102411138A (zh) 一种机器人声源定位方法
Ginn et al. Noise source identification techniques: simple to advanced applications
CN103278801A (zh) 一种变电站噪声成像侦测装置及侦测计算方法
RU2005105753A (ru) Система и способ для формирования луча с использованием микрофонной решетки
CN102147458B (zh) 一种针对宽带声源的波达方向估计方法及其装置
KR20160086831A (ko) 음장 재현 장치 및 방법, 그리고 프로그램
WO2013072554A2 (en) Spatial visual effect creation and display such as for a screensaver
WO2018003158A1 (ja) 相関関数生成装置、相関関数生成方法、相関関数生成プログラムおよび波源方向推定装置
JP2018091647A (ja) 信号処理装置、方位算出方法及び方位算出プログラム
Hosseini et al. Time difference of arrival estimation of sound source using cross correlation and modified maximum likelihood weighting function
Zhong et al. Direction of arrival tracking of an underwater acoustic source using particle filtering: Real data experiments
Lincke et al. Synthesizing coherence loss by atmospheric turbulence in virtual microphone array signals
Tourbabin et al. Speaker localization by humanoid robots in reverberant environments
JP5071938B2 (ja) 音源探査装置
CN102707261A (zh) 一种麦克风阵列声源定位系统
Buerger et al. The spatial coherence of noise fields evoked by continuous source distributions
Mansour et al. Blind separation of underwater acoustic signals

Legal Events

Date Code Title Description
ENP Entry into the national phase

Ref document number: 2018523292

Country of ref document: JP

Kind code of ref document: A

121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 17812908

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 17812908

Country of ref document: EP

Kind code of ref document: A1