WO2021079513A1 - 信号処理器、レーザレーダ装置及び風力タービン - Google Patents

信号処理器、レーザレーダ装置及び風力タービン Download PDF

Info

Publication number
WO2021079513A1
WO2021079513A1 PCT/JP2019/042018 JP2019042018W WO2021079513A1 WO 2021079513 A1 WO2021079513 A1 WO 2021079513A1 JP 2019042018 W JP2019042018 W JP 2019042018W WO 2021079513 A1 WO2021079513 A1 WO 2021079513A1
Authority
WO
WIPO (PCT)
Prior art keywords
wind speed
signal processor
unit
wind
value
Prior art date
Application number
PCT/JP2019/042018
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 PCT/JP2019/042018 priority Critical patent/WO2021079513A1/ja
Priority to JP2021553876A priority patent/JP6995256B2/ja
Publication of WO2021079513A1 publication Critical patent/WO2021079513A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P5/00—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft
    • G01P5/18—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the time taken to traverse a fixed distance
    • G01P5/20—Measuring speed of fluids, e.g. of air stream; Measuring speed of bodies relative to fluids, e.g. of ship, of aircraft by measuring the time taken to traverse a fixed distance using particles entrained by a fluid stream
    • 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
    • G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
    • G01S17/88—Lidar systems specially adapted for specific applications
    • G01S17/95—Lidar systems specially adapted for specific applications for meteorological use
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Definitions

  • the signal processor sends out a plurality of modulated light signals generated by using a laser light source in each of a plurality of directions in the observation space, and the plurality of modulated light signals generated by backscattering the plurality of modulated light signals.
  • a signal processor that operates in cooperation with a laser radar sensor that receives the scattered light signals of the above as a plurality of received light signals, detects the plurality of received light signals, and generates a plurality of received signals corresponding to the plurality of directions.
  • the wind condition measuring unit calculates the measured value of the line-of-sight velocity of the wind existing in each of the plurality of distance ranges from the laser radar sensor in the plurality of directions by frequency analysis of the plurality of received signals.
  • Data assimilation unit 74 the predicted wind field at the present time t k which is calculated previously, is assimilated to the estimated wind field at the present time t k which is supplied from the wind speed estimator 60 calculates the filtered wind field by.
  • This filter wave wind speed field is a one-dimensional discrete distribution consisting of the filter wave wind speed values V fil (k) (1), ..., V fil (k) (Nbin). That is, the data assimilation unit 74, the predicted wind speed value V est (k): (i) by assimilation to the estimated wind speed value V cen the present time t k (k) (i) filtered wind speed V fil the current time t k (k) Calculate (i).
  • the time evolution calculation unit 78 uses the time evolution equation for one-dimensional space derived from the Navier-Stokes equation to sequentially calculate the time evolution of the wind speed field from the filtered wind speed field, thereby causing the future time tk + 1 , ... , T k + C, the predicted wind speed field can be calculated. That is, the time evolution calculation unit 78 calculates the predicted wind speed values Vest (k + 1) (i), ..., Vest (k + C) (i) using the time evolution equation, and the predicted wind speed value calculated after the C cycle.
  • the vest (k + C) (i) can be output to the outside as the predicted wind speed data PWD. At this time, the predicted wind speed value Vest (k + 1) (i) calculated after one cycle is output to the delay element 79.
  • u is a wind speed vector
  • ⁇ is a viscosity coefficient
  • ⁇ is a density
  • p is a pressure
  • F is an external force vector
  • ⁇ is a time variable.
  • the wind speed vector u (n) as shown in the following equations (18.1) and (18.2) is used. Introduce an intermediate wind speed vector u *.
  • the superscript n is an integer representing the local time tau n between future time t k + C after the current time t k and C cycle.
  • the local time ⁇ n changes in units of time step d ⁇ .
  • the time evolution calculation unit 78 can output the calculated wind speed value u i ( ⁇ ) as the predicted wind speed value West (k + C) (i) at the future time tk + C.
  • the time evolution calculation unit 78 calculates a one cycle after the wind speed value u i ( ⁇ 0) by executing sequentially the iterative calculation using the above time evolution equation.
  • ⁇ 0 is the number of iterative calculations over one cycle, and is a value corresponding to the time closest to the time tk + 1 of the next cycle.
  • Time evolution calculation unit 78, the calculated wind velocity values u i ( ⁇ 0) can be output as the time for the next cycle t k + 1 of the predicted wind speed V est (k + 1) ( i).
  • the control unit 49 the cycle number k indicating the time t k is set to "1" to start the control of the wind condition measurement (step ST11).
  • the wind speed calculation unit 62 includes measured values Vobs (k) (i, 1) to Vobs (k) (i, J) and evaluation values q (k) (i, 1) to q. (k)
  • the estimated wind speed values V fen (k) (i) are calculated using (i, J) (step ST26).
  • step ST27 when the range bin number i has not reached the upper limit value Nbin (NO in step ST27), the effectiveness evaluation unit 61 increments the range bin number i by 1 (step ST28) and performs steps ST23 to ST25. Execute. Further, the wind speed calculation unit 62 executes step ST26. Finally, when the number i reaches the upper limit value Nbin (YES in step ST27), the estimated wind speed values V fen (k) (i) are calculated for all the values of the range bin number i.
  • the data assimilation unit 74 first sets the range bin number i to “-1” (step ST). Next, the data assimilation unit 74 executes steps ST43 to ST49 when the cycle number k is “1” or less (NO in step ST42), and when the cycle number k is larger than “1” (step ST42). YES of ST42), steps ST51 to ST56 are executed.
  • step ST42 the cycle number k is “1” or less (NO in step ST42) will be described.
  • step ST43 data assimilation unit 74 sets filtered wind speed V fil a (k) (i) to "0" (step ST44) .
  • step ST45 data assimilation unit 74 sets the filtered wind speed V fil (k) (i) the estimated wind speed value V cen (k) (i) ( step ST46).
  • i up is the first valid estimated wind speed value V fen (k) (i up ) found when searching for the estimated wind speed values V fen (k) (i) in ascending order of range bin number i. ),
  • i dn is the first valid estimated wind speed value V fen (k) (i ) found when searching for the estimated wind speed values V fen (k) (i) in descending order of the range bin number i. dn ) range bin number.
  • step ST48 the data assimilation unit 74 increments the range bin number i by 1 if the range bin number i has not reached the upper limit value Nbin (NO in step ST48). (Step ST49), and the processing after step ST43 is executed. Finally, when the range bin number i reaches the upper limit value Nbin (YES in step ST48), the data assimilation process is completed.
  • the signal processor 41 of the first embodiment is distributed in the central axis direction of the laser radar sensor 11 based on the measured values of the line-of-sight velocity of the wind in the plurality of line-of-sight directions LOS 1 to LOS J1.
  • the estimated wind velocity field of the dimension is calculated, and the time evolution of the wind velocity field is calculated from the estimated wind velocity field of the one dimension using the time evolution equation for the one-dimensional space derived from the Navier-Stokes equation. This makes it possible to calculate the predicted wind speed field with a small amount of calculation while maintaining high accuracy.
  • All or part of the functions of the signal processor 41 described above are, for example, a semiconductor integrated circuit having a DSP (Digital Signal Processor), an ASIC (Application Specific Integrated Circuit), or an FPGA (Field-Programmable Gate Array). It can be realized by the processor of. Alternatively, all or part of the functions of the signal processor 41 may be performed by one or more processors including an arithmetic unit such as a CPU (Central Processing Unit) or a GPU (Graphics Processing Unit) that executes software or firmware program code. It may be realized.
  • arithmetic unit such as a CPU (Central Processing Unit) or a GPU (Graphics Processing Unit) that executes software or firmware program code. It may be realized.
  • FIG. 11 is a block diagram showing a schematic configuration of the wind condition estimation unit 47 according to the present embodiment.
  • the wind condition estimation unit 48 has a wind speed estimation unit 60 and a wind speed prediction unit 73.
  • the wind speed prediction unit 73 includes a data assimilation unit 76, a time step setting unit 77, a time evolution calculation unit 78, and a delay element 79.
  • the configuration of the wind condition estimation unit 48 of the present embodiment is the same as the configuration of the wind condition estimation unit 46 of the first embodiment except that the data assimilation unit 76 is provided instead of the data assimilation unit 74 of FIG. is there.
  • the data assimilation unit 76 is supplied with the number of valid measured values N q (k) (i) from the wind speed estimation unit 60.
  • Data assimilation unit 76 based on the number of valid measurements N q (k) (i) , should be assigned to the estimated wind field weighting factor ⁇ and the current time t k to be assigned to the prediction wind field at the current time t k It is characterized in that the weighting coefficient 1- ⁇ is adaptively set.
  • the laser radar device 4 includes a laser radar sensor 12 having a configuration different from that of the laser radar sensor 11 and a signal processor 41.
  • the configuration of the signal processor 41 of the present embodiment is the same as the configuration of the signal processor 41 of the first embodiment.
  • the signal processor 41 of the present embodiment either the signal processor 42 of the second embodiment or the signal processor 43 of the third embodiment may be used.
  • the tunable laser light source 22 of the present embodiment has a function of switching the wavelengths of the output laser light ⁇ 1 , ⁇ 2 , ..., ⁇ J according to the control signal supplied from the signal processor 41. There is.
  • the optical transceiver 32 includes a wavelength separation device 34 such as an optical multiplexer is configured to include an optical antenna 35 1 ⁇ 35 J.
  • Wavelength demultiplexing element 34 depending on the wavelength of the input amplified light pulse, connects the one of the optical input-output end and the optical antenna 35 1 ⁇ 35 J of the optical circulator 26 is optically.
  • the wavelength separating element 34 outputs the amplified light pulse only to the j-th optical antenna 35 j .
  • the amplified optical pulse is transmitted from the j-th optical antenna 35 j in the line-of-sight direction LOS j.
  • 1,4 laser radar device 11,12 laser radar sensor, 21 laser light source, 22 wavelength variable laser light source 22, 23 optical distributor, 24 optical modulator, 25 optical amplifier, 26 optical circulator, 27 optical combiner, 28 optical detector, 31 and 32 optical transceiver, 33 an optical switch, 34 wavelength separation device, 35 1 ⁇ 35 J optical antenna, 41, 42, and 43 signal processor, 44 wind measuring unit, 45 a buffer memory, 46 and 47, 48 Wind condition estimation unit, 49 Control unit, 53 Quality evaluation unit, 56 Spectrum integration unit, 57 Line-of-sight speed calculation unit, 60 Wind speed estimation unit, 61 Effectiveness evaluation unit, 62 Wind speed calculation unit, 71, 72, 73 Wind speed prediction unit , 74, 75, 76 Data assimilation unit, 77-hour step setting unit, 78-hour development calculation unit, 79 delay element, 90 signal processing circuit, 91 processor, 92 memory, 93 storage device, 94 input / output interface unit, 95 signal path , 100 wind turbine, 101 nacelle, 102 tower, 103 hub, 104A

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Electromagnetism (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Optical Radar Systems And Details Thereof (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

観測空間の複数方位にそれぞれ対応する複数の受信信号を生成するレーザレーダセンサ（１１）と連携して動作する信号処理器（４１）が開示される。信号処理器（４１）は、複数の受信信号を周波数解析することにより、レーザレーダセンサ（１１）からの複数の距離範囲にそれぞれ存在する風の視線方向速度の計測値を算出し、当該視線方向速度の計測値に基づいて、レーザレーダセンサ（１１）の中心軸方向に分布する１次元風速場を現時刻の推定風速場として算出する。信号処理器（４１）は、さらに、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、当該推定風速場からの風速場の時間発展を計算することにより将来時刻の予測風速場を算出する。

Description

信号処理器、レーザレーダ装置及び風力タービン
　本発明は、レーザレーダ技術に関し、特に、観測空間にて反射または散乱されたレーザ光に基づいて風速及び風向などの風況特性を計測するためのレーザレーダ技術に関する。
　レーザレーダ技術は、観測空間に照射された後に反射または散乱されたレーザ光に基づいて当該観測空間の性質及び状態を観測するための技術であり、ライダ（ＬＩｇｈｔ　Ｄｅｔｅｃｔｉｏｎ　Ａｎｄ　Ｒａｎｇｉｎｇ，ＬＩＤＡＲ）と呼ばれることもある。このようなレーザレーダ技術は、観測空間内に浮遊するエアロゾル（ａｅｒｏｓｏｌｓ）及び分子などの微粒子にて後方散乱されたレーザ光に基づいて風速及び風向などの風況特性（ｗｉｎｄ　ｐｒｏｐｅｒｔｉｅｓ）を計測するために広く使用されている。たとえば、風力発電の用途では、発電効率の向上のために、レーザレーダによりリアルタイムに計測された風況特性を用いて風力タービン（風車）のナセルの向きあるいはブレードのピッチ角を適応的に制御することが行われている。また風力発電システムまたはウィンドファーム（ｗｉｎｄ　ｆａｒｍ）の建設予定地の風況観測のために、レーザレーダ技術が活用されることもある。ここで、ウィンドファームとは、陸上または洋上に配置された複数基の風力発電システムを有する発電所をいう。
　レーザレーダを用いた風況特性の推定技術は、たとえば、下記の特許文献１及び非特許文献１に開示されている。特許文献１には、あらかじめ決められた気流モデルを用いて、ＬＩＤＡＲセンサの検出出力から、風力タービンのロータに到来する将来の風況場（ｗｉｎｄ　ｆｉｅｌｄ）の風況特性を予測する技術がいくつか開示されている。特許文献１に開示されている技術の１つは、テイラーの凍結乱流仮説（Ｔａｙｌｏｒ’ｓ　ｆｒｏｚｅｎ　ｔｕｒｂｕｌｅｎｃｅ　ｈｙｐｏｔｈｅｓｉｓ）に基づく気流モデルを用いて風況特性を予測する方法である。しかしながら、この方法は、少ない演算量で将来の風況特性を予測することができるものの、予測精度が低いという課題がある。
欧州特許第２５８１７６１号明細書（たとえば、段落［００５２］～［００５５］，［００５９］～［００６３］）
P. Towers and B. L. Jones, "Real-time wind field reconstruction from LiDAR measurements using a dynamic wind model and state estimation," Wind Energy, vol. 19, issue 1, pp. 133-150, 2016.
　風況特性の高精度予測技術としては、３次元空間及び時間に関するナビエ・ストークス方程式の直接数値シミュレーション（Ｄｉｒｅｃｔ　Ｎｕｍｅｒｉｃａｌ　Ｓｉｍｕｌａｔｉｏｎ，ＤＮＳ）が知られている。しかしながら、そのような直接数値シミュレーションは膨大な演算量を必要とし、演算負荷が極めて高いという課題がある。
　そこで非特許文献１には、ナビエ・ストークス方程式を簡易化することにより導出された風況推定モデルを用いて、風況特性の時間発展の数値計算を行う方法が提案されている。この方法は、前述の直接数値シミュレーション（ＤＮＳ）と比べれば、少ない演算量で数値計算を行うことが可能である。しかしながら、ハードウェアに実装する観点からは、さらに少ない演算量で数値計算を行うことができる方法が求められている。
　上記に鑑みて本発明の目的は、高い精度を保ちつつ少ない演算量で風況特性を予測することができる信号処理器、レーザレーダ装置及びこれを有する風力タービンを提供することである。
　本発明の一態様による信号処理器は、レーザ光源を用いて生成された複数の変調光信号を観測空間の複数方位にそれぞれ送出し、前記複数の変調光信号が後方散乱することにより生じた複数の散乱光信号を複数の受信光信号として受信し、前記複数の受信光信号を検波して前記複数方位にそれぞれ対応する複数の受信信号を生成するレーザレーダセンサと連携して動作する信号処理器であって、前記複数の受信信号を周波数解析することにより、前記複数方位において前記レーザレーダセンサからの複数の距離範囲にそれぞれ存在する風の視線方向速度の計測値を算出する風況計測部と、前記視線方向速度の計測値に基づいて、前記レーザレーダセンサの中心軸方向に分布する１次元風速場を現時刻の推定風速場として算出する風速推定部と、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、前記推定風速場からの風速場の時間発展を計算することにより将来時刻の予測風速場を算出する風速予測部とを備えることを特徴とする信号処理器。
　本発明の一態様によれば、複数方位における風の視線方向速度の計測値に基づいてレーザレーダセンサの中心軸方向に分布する１次元の推定風速場が算出され、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、当該１次元の推定風速場から風速場の時間発展が計算される。これにより、高い精度を保ちつつ少ない演算量で予測風速場を算出することが可能となる。
本発明に係る実施の形態１のレーザレーダ装置の概略構成を示すブロック図である。 実施の形態１の風況計測部の概略構成を示すブロック図である。 図３Ａ～図３Ｃは、各種信号波形を表すグラフである。 実施の形態１における風況推定部の概略構成を示すブロック図である。 レーザレーダ装置を有する風力タービンを例示する図である。 現時刻の予測風速値の分布と将来時刻の予測風速値の分布とを概念的に表すグラフである。 実施の形態１の信号処理器の動作手順の一例を概略的に示すフローチャートである。 実施の形態１に係るデータ同化処理の手順の一例を概略的に示すフローチャートである。 実施の形態１の信号処理器のハードウェア構成例を示すブロック図である。 本発明に係る実施の形態２のレーザレーダ装置における信号処理器の概略構成を示すブロック図である。 実施の形態２における風況推定部の概略構成を示すブロック図である。 実施の形態２に係るデータ同化処理の手順の一例を概略的に示すフローチャートである。 本発明に係る実施の形態３のレーザレーダ装置における信号処理器の概略構成を示すブロック図である。 実施の形態３における風況推定部の概略構成を示すブロック図である。 実施の形態３に係るデータ同化処理の手順の一例を概略的に示すフローチャートである。 本発明に係る実施の形態４のレーザレーダ装置の概略構成を示すブロック図である。
　以下、図面を参照しつつ、本発明に係る種々の実施の形態及びこれらの変形例について詳細に説明する。なお、図面全体において同一符号を付された構成要素は、同一構成及び同一機能を有するものとする。
実施の形態１．
　図１は、本発明に係る実施の形態１のレーザレーダ装置１の概略構成を示すブロック図である。図１に示されるようにレーザレーダ装置１は、レーザレーダセンサ１１と信号処理器４１とを備えて構成されたライダ（Ｌｉｇｈｔ　Ｄｅｔｅｃｔｉｏｎ　ａｎｄ　Ｒａｎｇｉｎｇ，ＬＩＤＡＲ）装置である。本実施の形態のレーザレーダ装置１は、風況特性としての風速を計測及び予測するように構成されている。
　レーザレーダセンサ１１は、レーザ光源を用いて生成された複数の変調光信号を、観測空間の複数方位すなわち視線方向ＬＯＳ１，ＬＯＳ２，…，ＬＯＳＪにそれぞれ送出し、当該複数の変調光信号が観測空間内に浮遊するエアロゾル（ａｅｒｏｓｏｌｓ）及び分子などの微粒子にて後方散乱することにより生じた複数の散乱光信号を複数の受信光信号として受信し、当該複数の受信光信号を検波して前記複数方位にそれぞれ対応する複数の受信信号を生成する機能を有している。ここで、Ｊは、あらかじめ定められた視線方向ＬＯＳ１～ＬＯＳＪの個数である。本実施の形態では、視線方向ＬＯＳ１～ＬＯＳＪの個数は３個以上であるが、これに限定されず、２個でもよい。
　具体的には、図１に示されるように、レーザレーダセンサ１１は、単一波長（単一周波数）の連続波であるレーザ光を出力するレーザ光源２１と、レーザ光源２１から入力されたレーザ光を送信光信号と局発光信号とに分配する光分配器２３と、光分配器２３から入力された送信光信号に周波数変調及び振幅変調を施すことにより一連の変調光パルス（複数の変調光信号）を生成する光変調器２４と、光変調器２４から入力された一連の変調光パルスを増幅して一連の増幅光パルスを生成する光増幅器２５と、光増幅器２５の光出力端と光学的に接続された光サーキュレータ２６と、光増幅器２５から光サーキュレータ２６を介して入力された一連の増幅光パルスを観測空間内の視線方向ＬＯＳ１，ＬＯＳ２，…，ＬＯＳＪに送出する光送受信部３１と、光合波器２７と、光検出器２８とを備えている。光合波器２７と光検出器２８とにより光ヘテロダイン検波器（ｏｐｔｉｃａｌ　ｈｅｔｅｒｏｄｙｎｅ　ｄｅｔｅｃｔｏｒ）が構成される。
　レーザ光源２１としては、半導体レーザ光源または固体レーザ光源が使用可能である。また光増幅器２５としては、音響光学変調器（Ａｃｏｕｓｔｏ－Ｏｐｔｉｃ　Ｍｏｄｕｌａｔｏｒ）が使用可能であるが、これに限定されるものではない。光増幅器２５は、長距離離れた遠方領域の風速を計測するために設けられており、たとえば、光ファイバ増幅器により構成可能である。
　光送受信部３１は、光サーキュレータ２６の光入出力端と光学的に結合された光スイッチ３３と、光アンテナ３５１～３５Ｊとを備えて構成されている。光スイッチ３３は、信号処理器４１の制御部４９から供給された制御信号に従って動作し、光サーキュレータ２６の光入出力端と光アンテナ３５１～３５Ｊのいずれかとの間を光学的に接続する。光スイッチ３３がｊ番目の光アンテナ３５ｊと光サーキュレータ２６の光入出力端との間を接続している間に、光サーキュレータ２６から光送受信部３１に入力された増幅光パルスは、ｊ番目の光アンテナ３５ｊから視線方向ＬＯＳｊに送出される。このとき、観測空間にて当該増幅光パルスが後方散乱することにより散乱光信号が発生すると、ｊ番目の光アンテナ３５ｊは当該散乱光信号を受信光信号として受信する。当該受信光信号は、光スイッチ３３を介して光サーキュレータ２６に伝搬する。
　光サーキュレータ２６は、非相反型の光学素子である。すなわち、光サーキュレータ２６は、光増幅器２５から順方向に入力された増幅光パルスを光送受信部３１に出力するが、当該増幅光パルスを光合波器２７に出力しない。一方、光サーキュレータ２６は、光送受信部３１から逆方向に入力された受信光信号を光合波器２７に出力するが、当該受信光信号を光増幅器２５に出力しない。このような光サーキュレータ２６は、たとえば、１／４波長板と偏光ビームスプリッタとにより構成可能である。
　光合波器２７は、光分配器２３から入力された局発光信号と、光送受信部３１から光サーキュレータ２６を介して入力された受信光信号とを合波して光ビート信号を生成し、当該光ビート信号を光検出器２８に出力する。光検出器２８は、光ビート信号を検波してアナログ受信信号（検波信号）Ｒｊ（ｈ，ｔ）を生成し、当該アナログ受信信号Ｒｊ（ｈ，ｔ）を信号処理器４１に出力する。ここで、アナログ受信信号Ｒｊ（ｈ，ｔ）は、ｊ番目の光アンテナ３５ｊを通じて受信されたｈ番目の受信光信号から得られた受信信号をいう。添え字ｊは、光アンテナ３５ｊの番号を示す１～Ｊの範囲内の整数であり、ｈは、受信光信号の番号（パルスヒット番号）を示す０～Ｈ－１の範囲内の整数であり（Ｈは正整数）、ｔは、時間を示す変数である。
　本実施の形態の光検出器２８は、バランスド受信器（ｂａｌａｎｃｅｄ　ｒｅｃｅｉｖｅｒ）構成を有しており、たとえば、フォトダイオードなどの一対の受光素子により構成可能である。
　レーザレーダセンサ１１の内部における光経路は、光ファイバにより構成可能である。すなわち、レーザ光源２１と光分配器２３との間の光経路、光分配器２３と光変調器２４との間の光経路、光変調器２４と光増幅器２５との間の光経路、光増幅器２５と光サーキュレータ２６との間の光経路、光サーキュレータ２６と光送受信部３１との間の光経路、光分配器２３と光合波器２７との間の光経路、光サーキュレータ２６と光合波器２７との間の光経路は、光ファイバにより構成可能である。ただし、これら光経路は、光ファイバに限定されるものではない。
　図１に示されるように信号処理器４１は、風況計測部４４、バッファメモリ４５、風況推定部４６及び制御部４９を備える。制御部４９は、レーザレーダセンサ１１におけるレーザ光源２１，光変調器２４，光増幅器２５及び光スイッチ３３の動作をそれぞれ制御するための制御信号を供給し、風況計測部４４，バッファメモリ４５及び風況推定部４６の動作をそれぞれ制御するための制御信号ＭＣ，ＢＣ，ＥＣを生成する。
　風況計測部４４は、アナログ受信信号Ｒｊ（ｈ，ｔ）を周波数解析することにより、視線方向ＬＯＳ１～ＬＯＳＪにおいてレーザレーダセンサ１１からの複数の距離範囲Δ（ｉ＝１），Δ（ｉ＝２），…，Δ（ｉ＝Ｎｂｉｎ）にそれぞれ存在する風（流体）の視線方向速度の計測値を算出する。ここで、ｉは、距離範囲Δ（ｉ）に割り当てられた番号（レンジビン番号）であり、Ｎｂｉｎは、３以上の整数である。
　図２は、実施の形態１の風況計測部４４の概略構成を示すブロック図である。図２に示されるように風況計測部４４は、時間ゲーティング部５１、領域変換部５２、品質評価部５３、スペクトル積算部５６及び視線速度算出部５７を有している。
　時間ゲーティング部５１は、各アナログ受信信号Ｒｊ（ｈ，ｔ）に対してＡ／Ｄ変換と時間ゲーティングとを実行することにより、距離範囲Δ（１），Δ（２），…，Δ（Ｎｂｉｎ）にそれぞれ対応する時間領域信号ｒ１，ｊ（ｈ，τ），ｒ２，ｊ（ｈ，τ），…，ｒＮｂｉｎ，ｊ（ｈ，τ）を生成する。ここで、τは、各距離範囲Δ（ｉ）に対応する時間領域内の局所的な時刻を表すサンプル番号である。
　たとえば、時間ゲーティング部５１は、各アナログ受信信号Ｒｊ（ｈ，ｔ）にＡ／Ｄ変換を施してディジタル信号を生成し、当該ディジタル信号を複数の時間領域信号ｒ１，ｊ（ｈ，τ）～ｒＮｂｉｎ，ｊ（ｈ，τ）に分割することができる。あるいは、時間ゲーティング部５１は、Ａ／Ｄ変換と時間ゲーティングとを同時に実行することにより、時間領域信号ｒ１，ｊ（ｈ，τ）～ｒＮｂｉｎ，ｊ（ｈ，τ）を生成してもよい。なお、レーザレーダセンサ１１における光検出器２８は、アナログ受信信号Ｒｊ（ｈ，ｔ）をディジタル受信信号に変換するＡ／Ｄ変換器を有していてもよい。この場合、時間ゲーティング部５１は、Ａ／Ｄ変換を行う機能を有していなくてもよい。
　図３Ａは、ｊ番目の光アンテナ３５ｊを通じて受信されたＨ個の受信光信号から得られたアナログ受信信号Ｒｊ（０，ｔ），Ｒｊ（１，ｔ），…，Ｒｊ（Ｈ－１，ｔ）の振幅波形を概念的に表すグラフである。図３Ａに示されるように、各アナログ受信信号Ｒｊ（ｈ，ｔ）から、時間幅Δｔの時間領域における信号成分が選択され、当該信号成分から時間領域信号ｒｉ,ｊ（ｈ，τ）が生成される。
　次に、図２を参照すると、領域変換部５２は、時間領域信号ｒｉ,ｊ（ｈ，τ）に対し、サンプル番号τについて時間領域から周波数領域への領域変換処理を実行することにより、周波数領域信号ｓｉ,ｊ（ｈ，ｆ）を生成する。ここで、ｆは、周波数を表すサンプル番号（周波数番号）である。領域変換処理としては、高速フーリエ変換（Ｆａｓｔ　Ｆｏｕｒｉｅｒ　Ｔｒａｎｓｆｏｒｍ，ＦＦＴ）などの離散直交変換が挙げられる。図３Ｂは、高速フーリエ変換（ＦＦＴ）により得られる周波数領域信号ｓｉ,ｊ（０，ｆ），ｓｉ,ｊ（１，ｆ），…，ｓｉ,ｊ（Ｈ－１，ｆ）のスペクトルを概念的に表すグラフである。
　周波数領域信号ｓｉ,ｊ（ｈ，ｆ）は、品質評価部５３とスペクトル積算部５６とに供給される。品質評価部５３は、周波数領域信号ｓｉ,ｊ（ｈ，ｆ）の信号品質を評価し、その評価結果である品質評価値Ｑｉ,ｊ（ｈ）をスペクトル積算部５６に出力する。
　具体的には、品質評価部５３は、周波数領域信号ｓｉ,ｊ（ｈ，ｆ）から信号対雑音比（ｓｉｇｎａｌ－ｔｏ－ｎｏｉｓｅ　ｒａｔｉｏ，ＳＮＲ）を算出し、ＳＮＲが閾値以上であれば、品質良好を示す値（たとえば「１」）を有する品質評価値Ｑｉ,ｊ（ｈ）をスペクトル積算部５６に出力し、ＳＮＲが閾値未満であれば、品質不良を示す値（たとえば「０」）を有する品質評価値Ｑｉ,ｊ（ｈ）をスペクトル積算部５６に出力することができる。
　スペクトル積算部５６は、品質評価値Ｑｉ,ｊ（ｈ）に基づき、品質良好と評価された周波数領域信号ｓｉ,ｊ（ｈ，ｆ）のスペクトルのみをパルスヒット番号ｈについて積算（インコヒーレント積分）することにより積算スペクトル｜ｓｉ,ｊ（ｆ）｜２を算出することができる。図３Ｃは、積算スペクトル｜ｓｉ,ｊ（ｆ）｜２を概念的に表すグラフである。視線速度算出部５７は、積算スペクトル｜ｓｉ,ｊ（ｆ）｜２に現れる極大ピークの周波数番号ｆＤを検出し、当該検出された周波数番号ｆＤからドップラ周波数を算出する。視線速度算出部５７は、当該算出されたドップラ周波数から、各視線方向ＬＯＳｊにおける各距離範囲Δ（ｉ）に存在する風の視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）を算出することができる。ここで、上付き添え字ｋは、計測時刻ｔｋを表す整数である。
　またスペクトル積算部５６は、積算スペクトルの算出に使用された周波数領域信号ｓｉ,ｊ（ｈ，ｆ）の個数から１を減算して得た値すなわち積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）を出力する。この積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）は、視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の品質を表す値である。積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）が一定値以上であれば、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の品質は良好であり、積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）が低すぎれば、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の品質は不良であると評価することができる。
　風況計測部４４は、視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）と積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）との組合せを風計測データＷＭＤとしてバッファメモリ４５に出力する。バッファメモリ４５は、風計測データＷＭＤを一時的に記憶する。
　次に、図１の風況推定部４６の構成について説明する。
　図４は、実施の形態１の風況推定部４６の概略構成を示すブロック図である。図４に示されるように風況推定部４６は、風速推定部６０及び風速予測部７１を有している。
　風速推定部６０は、バッファメモリ４５から風計測データＷＭＤを読み出し、風計測データＷＭＤに含まれる視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）及び積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）に基づいて、レーザレーダセンサ１１の中心軸方向に分布する１次元風速場を現時刻ｔｋの推定風速場として算出することができる。現時刻ｔｋの推定風速場は、推定風速値Ｖｃｅｎ (ｋ)（１），…，Ｖｃｅｎ (ｋ)（Ｎｂｉｎ）からなる１次元の離散的な分布である。具体的には、風速推定部６０は、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，Ｊ）を用いた非線形最小自乗法により推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出することができる。
　図５は、レーザレーダ装置１を有する風力タービン１００を例示する図である。図５に示される風力タービン１００は、レーザレーダ装置１を搭載するナセル１０１と、ナセル１０１を下方から支持するタワー１０２と、ナセル１０１に取り付けられたハブ１０３と、ハブ１０３を介してナセル１０１内のロータ軸（図示せず）と回転自在に連結されたブレード１０４Ａ，１０４Ｂ，１０４Ｃからなるロータとを備えて構成されている。レーザレーダ装置１は、ロータの回転面を介して観測空間に増幅光パルスを送出することができるようにナセル１０１上に配置されている。ナセル１０１は、風向きに追従させるためにタワー１０２の向きΨを変化させるヨー制御と、到来風速に応じた受風量の調整のためにブレード１０４Ａ～１０４Ｃの取付角Φを変化させるピッチ制御とを行うことができる制御装置（図示せず）を搭載している。
　なお、図５には、ロータの回転面が風上側に位置するアップウィンド方式の風力タービン１００が示されているが、これに限定されるものではない。ロータの回転面が風下側に位置するダウンウィンド方式の風力タービンがレーザレーダ装置１を搭載する実施の形態もありうる。
　図５の例の場合、レーザレーダ装置１のレーザレーダセンサ１１の中心軸ＣＡの周りに４個の視線方向ＬＯＳ１～ＬＯＳ４が設けられている。レーザレーダ装置１は、視線方向ＬＯＳ１～ＬＯＳ４に増幅光パルス（増幅された変調光パルス）を送出し、各距離範囲Δ（ｉ）にて後方散乱された散乱光に基づいて視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，４）を算出することができる。図４の風速推定部６０は、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，４）及び積算回数Ｎｓｕｍ (ｋ)（ｉ，１）～Ｎｓｕｍ (ｋ)（ｉ，４）に基づいて、距離範囲Δ（ｉ）における推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出することができる。
　風速推定部６０は、有効性評価部６１及び風速算出部６２を有する。有効性評価部６１は、風計測データＷＭＤに含まれる積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）を用いて、視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の有効性を評価し、その評価結果である評価値ｑ(ｋ)（ｉ，ｊ）を風速算出部６２に出力する。本実施の形態の評価値ｑ(ｋ)（ｉ，ｊ）は、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）が有効であることを示す値（たとえば「１」）、もしくは、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）が無効であることを示す値（たとえば「０」）のいずれか一方である。
　たとえば、次式（１）に示されるように、有効性評価部６１は、積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）が閾値Ｎｔｈ以上であれば、評価値ｑ(ｋ)（ｉ，ｊ）を「１」に設定し、積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）が閾値Ｎｔｈ未満であれば、評価値ｑ(ｋ)（ｉ，ｊ）を「０」に設定することができる。

Figure JPOXMLDOC01-appb-I000001
　なお、本実施の形態では、有効性評価部６１は、積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）に基づいて評価値ｑ(ｋ)（ｉ，ｊ）を設定しているが、この代わりに、品質評価値Ｑｉ,ｊ（ｈ）に基づいて評価値ｑ(ｋ)（ｉ，ｊ）を設定してもよい。この場合には、図２の風況計測部４４は、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）と積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）との組合せに代えて、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）と品質評価値Ｑｉ,ｊ（ｈ）（ｈ＝０～Ｈ－１）との組合せを風計測データＷＭＤとして出力すればよい。
　風速算出部６２は、評価値ｑ(ｋ)（ｉ，ｊ）に基づき、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）のうち各距離範囲Δ（ｉ）に対応する有効な計測値の個数Ｎｑ (ｋ)（ｉ）が閾値以上であるか否かを判定する。たとえば、個数Ｎｑ (ｋ)（ｉ）は、次式（２）により算出可能である。

Figure JPOXMLDOC01-appb-I000002
　風速算出部６２は、有効な計測値の個数Ｎｑ (ｋ)（ｉ）が閾値以上であると判定したときは、当該有効な計測値に基づいて推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出する。有効な計測値の個数Ｎｑ (ｋ)（ｉ）が閾値未満であると判定したときは、風速算出部６２は、推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出しない。これにより、精度の低い推定風速値Ｖｃｅｎ (ｋ)（ｉ）の算出を回避することができる。
　風速算出部６２は、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，Ｊ）及び評価値ｑ(ｋ)（ｉ，１）～ｑ(ｋ)（ｉ，Ｊ）を用いた非線形最小自乗法により推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出することができる。この非線形最小自乗法の一例を以下に説明する。説明の便宜上、４つの視線方向ＬＯＳ１～ＬＯＳ４の場合について説明するが、これに限定されるものではない。
　先ず、図５に示されるように、視線方向ＬＯＳ１～ＬＯＳ４における距離範囲Δ（ｉ）に存在する風の速度ベクトルＶｐ (k)（ｉ，１），Ｖｐ (k)（ｉ，２），Ｖｐ (k)（ｉ，３），Ｖｐ (k)（ｉ，４）の方向が中心軸ＣＡの方向と平行であると仮定する。この仮定の下では、速度ベクトルＶｐ (k)（ｉ，１），Ｖｐ (k)（ｉ，２），Ｖｐ (k)（ｉ，３），Ｖｐ (k)（ｉ，４）の大きさ（絶対値）は、次式（３．１）～（３．４）で与えられる。

Figure JPOXMLDOC01-appb-I000003
　式（３．１）～（３．４）において、ｖｃｅｎ (ｋ)（ｉ）は、中心軸ＣＡにおけるｉ番目の距離範囲Δ（ｉ）の風速値であり、δｖ（ｉ）は、垂直シア（鉛直方向における２点間の風速変化率）であり、δｈ（ｉ）は、水平シア（水平方向における２点間の風速変化率）である。Ｌ（ｉ）は次式（４）で与えられる。

Figure JPOXMLDOC01-appb-I000004
　式（４）において、Ｄ（ｉ）は、レーザレーダセンサ１１から距離範囲Δ（ｉ）までの距離であり、θは、中心軸ＣＡと視線方向ＬＯＳｊとがなす角度である。
　また、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，Ｊ）と速度ベクトルＶｐ (k)（ｉ，１）～Ｖｐ (k)（ｉ，４）との間の関係について次式（５．１）～（５．４）が成立する。

Figure JPOXMLDOC01-appb-I000005

　式（５．１）～（５．４）において、ｅ（ｉ，１）～ｅ（ｉ，４）は、誤差である。
　次に、次式（６．１）～（６．３）に示すようなベクトルＶｏｂｓ，ｖｃｅｎ，ｅを定義する（ここで、上付き添え字「Ｔ」は転置を示す）。

Figure JPOXMLDOC01-appb-I000006
　さらに、次式（７），（８）に示すような行列Ｑ，Ｌを定義する。

Figure JPOXMLDOC01-appb-I000007

Figure JPOXMLDOC01-appb-I000008
　ベクトルＶｏｂｓ，ｖｃｅｎ，ｅ及び行列Ｑ，Ｌを用いれば、式（５．１）～（５．４）は、次式（９）で表される。

Figure JPOXMLDOC01-appb-I000009
　非線形最小自乗法によれば、次式（１０）で表される誤差の二乗和ＳＳを最小化するベクトルｖｃｅｎの推定量を求めることが可能である。

Figure JPOXMLDOC01-appb-I000010
　非線形最小自乗法によれば、ベクトルｖｃｅｎの推定量Ｖｃｅｎは、次式（１１）で与えられる。

Figure JPOXMLDOC01-appb-I000011
　式（１１）において、Ｖｃｅｎ (ｋ)（ｉ）はｉ番目の距離範囲Δ（ｉ）における推定風速値、δＨ（ｉ）は水平シア推定値、δＶ（ｉ）は垂直シア推定値である。
　したがって、風速算出部６２は、式（１１）に従い、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，Ｊ）及び評価値ｑ(ｋ)（ｉ，１）～ｑ(ｋ)（ｉ，Ｊ）を用いて推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出することができる。計測環境が理想的である場合には、風速算出部６２は、評価値ｑ(ｋ)（ｉ，１）～ｑ(ｋ)（ｉ，Ｊ）を用いずに推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出することが可能である。
　次に、図４を参照しつつ風速予測部７１について説明する。風速予測部７１は、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、現時刻ｔｋの推定風速場からの風速場の時間発展を計算することにより将来時刻ｔｋ＋Ｃの予測風速場を算出することができる。ここで、Ｃは、サイクル回数を示す１以上の正整数であり、任意の正整数値に設定可能である。後述するように、１サイクル毎に、時間発展方程式を用いた反復計算が複数回実行される。現時刻ｔｋの推定風速場は、推定風速値Ｖｃｅｎ (ｋ)（１），…，Ｖｃｅｎ (ｋ)（Ｎｂｉｎ）からなる１次元の離散的な分布であり、将来時刻ｔｋ＋Ｃの予測風速場は、予測風速値Ｖｅｓｔ (ｋ＋Ｃ)（１），…，Ｖｅｓｔ (ｋ＋Ｃ)（Ｎｂｉｎ）からなる１次元の離散的な分布である。
　図４に示されるように風速予測部７１は、データ同化部７４、時間ステップ設定部７７、時間発展計算部７８及び遅延素子７９を有する。遅延素子７９は、以前に算出された現時刻ｔｋの予測風速値Ｖｅｓｔ (ｋ)（ｉ）（ｉ＝１～Ｎｂｉｎ）をデータ同化部７４に供給する。
　データ同化部７４は、以前に算出された現時刻ｔｋの予測風速場を、風速推定部６０から供給された現時刻ｔｋの推定風速場に同化させて濾波風速場を算出する。この濾波風速場は、濾波風速値Ｖｆｉｌ (ｋ)（１），…，Ｖｆｉｌ (ｋ)（Ｎｂｉｎ）からなる１次元の離散的な分布である。すなわち、データ同化部７４は、現時刻ｔｋの予測風速値Ｖｅｓｔ (ｋ)（ｉ）を現時刻ｔｋの推定風速値Ｖｃｅｎ (ｋ)（ｉ）に同化させて濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を算出する。
　たとえば、データ同化部７４は、現時刻ｔｋの予測風速場と現時刻ｔｋの推定風速場との重み付け和を濾波風速場として生成することができる。すなわち、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）は、重み係数α，１－αを用いて、次式（１２）で表される重み付け和として算出可能である。

Figure JPOXMLDOC01-appb-I000012
　図１の制御部４９は、重み係数α，１－αを計測環境に適した値に設定することができる。
　時間ステップ設定部７７は、濾波風速場に基づいて、時間発展計算部７８にて使用される時間ステップｄτを適応的に設定する。具体的には、時間ステップ設定部７７は、濾波風速場を構成する濾波風速値Ｖｆｉｌ (ｋ)（１）～Ｖｆｉｌ (ｋ)（Ｎｂｉｎ）のうちの最大風速値に基づいて、後述のクーラン条件（Ｃｏｕｒａｎｔ－Ｆｒｉｅｄｒｉｃｈｓ－Ｌｅｗｙ　Ｃｏｎｄｉｔｉｏｎ）を満たすように時間ステップｄτを適応的に設定することが好ましい。クーラン条件は、時間発展方程式の安定した数値解を確保するための条件である。この点については後述する。
　時間発展計算部７８は、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、濾波風速場からの風速場の時間発展を逐次的に計算することにより将来時刻ｔｋ＋１，…，ｔｋ＋Ｃにおける予測風速場を算出することができる。すなわち、時間発展計算部７８は、時間発展方程式を用いて予測風速値Ｖｅｓｔ (ｋ＋１)（ｉ），…，Ｖｅｓｔ (ｋ＋Ｃ)（ｉ）を算出し、Ｃサイクル後に算出された予測風速値Ｖｅｓｔ (ｋ＋Ｃ)（ｉ）を予測風速データＰＷＤとして外部に出力することができる。このとき、１サイクル後に算出された予測風速値Ｖｅｓｔ (ｋ＋１)（ｉ）は、遅延素子７９に出力される。図６は、現時刻ｔｋの予測風速値Ｖｅｓｔ (ｋ)（ｉ）の分布と、現時刻ｔｋから１サイクル後の将来時刻ｔｋ＋１における予測風速値Ｖｅｓｔ (ｋ＋１)（ｉ）の分布とを概念的に表すグラフである。ｉ＝－１の場合の予測風速値Ｖｅｓｔ (ｋ)（ｉ），Ｖｅｓｔ (ｋ＋１)（ｉ）は、レーザレーダセンサ１１に到来する風速値を表している。
　次に、時間発展方程式の導出法について説明する。
　通常速度の風は非圧縮流体とみなされる。非圧縮流体に関するナビエ・ストークス方程式は、次式（１３．１），（１３．２）で表現される。

Figure JPOXMLDOC01-appb-I000013
　ここで、ｕは風速ベクトル、νは粘性係数、ρは密度、ｐは圧力、Ｆは外力ベクトル、τは時間変数である。
　粘性係数νは非常に小さいので、風の粘性を考慮しなくてよい。外力ベクトルＦの項も無視することができる。よって、風の粘性が考慮されず、かつ外力が存在しないとの条件の下では、式（１３．１），（１３．２）から次式（１４．１），（１４．２）が導出される。

Figure JPOXMLDOC01-appb-I000014
　今、観測空間における風の流れが中心軸ＣＡに関して対称であると仮定する。また、中心軸方向における距離をｘとし、中心軸方向とは垂直な径方向における距離をｙとする円筒座標系（ｘ，ｙ）を考える。このとき、風速ベクトルｕは、次式（１５）で与えられる。

Figure JPOXMLDOC01-appb-I000015
　風速ベクトルｕのｘ成分ｕｘ及びｙ成分ｕｙを用いれば、式（１４．１），（１４．２）は、次式（１６．１），（１６．２），（１６．３）で表現される。

Figure JPOXMLDOC01-appb-I000016
　式（１６．３）を用いれば、式（１６．１），（１６．２）から次式（１７．１），（１７．２）が導出可能である。

Figure JPOXMLDOC01-appb-I000017
　次に、式（１７．１），（１７．２）を時間方向に離散化するために、次式（１８．１），（１８．２）に示されるような風速ベクトルｕ(ｎ)と中間的な風速ベクトルｕ＊とを導入する。

Figure JPOXMLDOC01-appb-I000018
　ここで、上付き添え字ｎは、現時刻ｔｋとＣサイクル後の将来時刻ｔｋ＋Ｃとの間の局所的な時刻τｎを表す整数である。局所的な時刻τｎは、時間ステップｄτの単位で変化する。
　風速ベクトルｕ(ｎ)を用いて式（１７．１），（１７．２）を時間方向に離散化すれば、次式（１９．１），（１９．２）が導出される。

Figure JPOXMLDOC01-appb-I000019
　中間的な風速ベクトルｕ＊を用いれば、式（１９．１）から次式（２０．１），（２０．２）が導出可能であり、式（１９．２）から次式（２１．１），（２１．２）が導出可能である。

Figure JPOXMLDOC01-appb-I000020

Figure JPOXMLDOC01-appb-I000021
　式（２０．２）の両辺を変数ｘで偏微分すれば、次式（２２．１）が導出される。また、式（２１．２）の両辺を変数ｙで偏微分すれば、次式（２２．２）が導出される。

Figure JPOXMLDOC01-appb-I000022
　連続の式（１６．３）を考慮すれば、式（２２．１），（２２．２）から次式（２３）が導出可能である。

Figure JPOXMLDOC01-appb-I000023
　観測空間の径方向における風速の変化率が極めて小さいと仮定すれば、変数ｙに関する１階偏微分の項を無視することができる。そこで、変数ｙに関する１階偏微分の項が考慮されないとの条件を採用すると、式（２０．１）から次の近似式（２４）が導出可能である。

Figure JPOXMLDOC01-appb-I000024
　また式（２０．２）を変形することで次式（２５）が導出される。

Figure JPOXMLDOC01-appb-I000025
　さらに、変数ｙに関する圧力ｐの２階偏微分が圧力ｐの線形関数であるとの条件を採用する。この条件の下でその線形関数をγｐ（γは定数）と表現すれば、式（２３）から次の近似式（２６）が導出可能である。

Figure JPOXMLDOC01-appb-I000026
　次に、式（２４），（２５），（２６）を空間方向に離散化するために、次の簡易表現式（２７．１），（２７．２），（２７．３）を導入する（ここで、ｄｘは、各距離範囲Δ（ｉ）の幅である）。

Figure JPOXMLDOC01-appb-I000027
　ところで、時間ステップｄτの大きさには制約が存在する。安定した数値解を得るためには、クーラン条件（Ｃｏｕｒａｎｔ－Ｆｒｉｅｄｒｉｃｈｓ－Ｌｅｗｙ　Ｃｏｎｄｉｔｉｏｎ）と呼ばれる次の条件式（２８）を満たすように時間ステップｄτを決定することが好ましい。

Figure JPOXMLDOC01-appb-I000028
　式（２８）において、βは１以上の余裕度係数である。条件式（２８）を常に満たすために、式（２８）の左辺の風速値ｕｉ (ｎ)は、風速値ｕ１ (ｎ)～ｕＮｂｉｎ (ｎ)のうちの最大風速値であることが好ましい。
　式（２４），（２６），（２５）を空間方向に離散化することで、それぞれ次式（２９），（３０），（３１）が導出される。

Figure JPOXMLDOC01-appb-I000029

Figure JPOXMLDOC01-appb-I000030

Figure JPOXMLDOC01-appb-I000031
　式（３０），（３１）において、Ｐｉ＝ｐｉ×ｄτ／ρ、である。
　空間方向の離散化の方法は、式（２９），（３０），（３１）の方法に限定されるものではない。安定した数値解を得るために、風上差分スキーム（Ｕｐｗｉｎｄ　Ｓｃｈｅｍｅ）と呼ばれる離散化方法が適用されてもよい。風上差分スキームは、空間方向の離散化の際に、波動進行方向反対側のサンプル差分を用いる離散化方法であり、β＝１の場合のクーラン条件が成立すれば、時間発展方程式の数値計算の安定化が保証される。式（２４）に風上差分スキームを適用すれば、次式（３２）が得られる。

Figure JPOXMLDOC01-appb-I000032
　式（３２）のＵｉ，Ｑｉは、次式（３３．１），（３３．２）により定義される。

Figure JPOXMLDOC01-appb-I000033
　式（３０）を用いてＰｉ＋１，Ｐｉを算出するために行列演算が適用できる。Ｐ０＝Ｐ１及びＰＫ＋１＝ＰＫと仮定すれば（Ｋは整数）、式（３０）は、次式（３４）の形に変形され得る。

Figure JPOXMLDOC01-appb-I000034
　式（３４）において、行列Ａは、次式（３５）で与えられる。

Figure JPOXMLDOC01-appb-I000035
　式（３５）において、Ｉ０は、すべての対角成分が１で、かつすべての非対角成分が零の行列である。
　したがって、次式（３６）に示されるように、行列Ａの逆行列Ａ－１を用いて、Ｐｉ＋１，Ｐｉを算出することができる。

Figure JPOXMLDOC01-appb-I000036
　以上により、時間発展計算部７８は、式（２９），（３６），（３１）の組合せ、または、式（３２），（３６），（３１）の組合せからなる時間発展方程式を用いて、１回の反復計算で、時刻τｎの風速値ｕｉ (ｎ)，ｕｉ＋１ (ｎ)から、時刻τｎ＋１（＝τｎ＋ｄτ）の風速値ｕｉ (ｎ＋１)を算出することができる。時間発展計算部７８は、１サイクル毎に、時間発展方程式を用いた反復計算を逐次的にκ回実行することにより、時刻τκ（＝τ０＋κ×ｄτ）における風速値ｕｉ (κ)を算出することが可能である。ここで、κは、反復計算の回数である。
　時間発展計算部７８は、データ同化部７４から濾波風速値Ｖｆｉｌ (ｋ)（ｉ）が与えられると、風速値ｕｉ (ｎ)の初期値ｕｉ (０)を濾波風速値Ｖｆｉｌ (ｋ)（ｉ）に設定する（ｉ＝１～Ｎｂｉｎ）。次に、時間発展計算部７８は、上記の時間発展方程式を用いた反復計算を逐次的に実行することによりＣサイクル後の風速値ｕｉ (Φ)を算出する。ここで、Φは、Ｃサイクルに亘る反復計算の回数であり、Ｃサイクル後の時刻ｔｋ＋Ｃに最も近い時刻に対応する値である。時間発展計算部７８は、算出された風速値ｕｉ (Φ)を、将来時刻ｔｋ＋Ｃの予測風速値Ｖｅｓｔ (ｋ＋Ｃ)（ｉ）として出力することができる。また時間発展計算部７８は、上記の時間発展方程式を用いた反復計算を逐次的に実行することにより１サイクル後の風速値ｕｉ (Φ０)を算出する。ここで、Φ０は、１サイクルに亘る反復計算の回数であり、次サイクルの時刻ｔｋ＋１に最も近い時刻に対応する値である。時間発展計算部７８は、算出された風速値ｕｉ (Φ０)を、次サイクルの時刻ｔｋ＋１の予測風速値Ｖｅｓｔ (ｋ＋１)（ｉ）として出力することができる。
　次に、図７のフローチャートを参照しつつ、上記した信号処理器４１の動作手順について説明する。図７は、信号処理器４１の動作手順の一例を概略的に示すフローチャートである。
　図７を参照すると、制御部４９は、時刻ｔｋを示すサイクル番号ｋを「１」に設定して風況計測の制御を開始する（ステップＳＴ１１）。
　次に、風況計測部４４は、レーザレーダセンサ１１から供給されたアナログ受信信号Ｒｊ（ｈ，ｔ）を基に風況計測を行う（ステップＳＴ１２）。すなわち、風況計測部４４は、上記のとおり、アナログ受信信号Ｒｊ（ｈ，ｔ）を周波数解析することにより、視線方向ＬＯＳ１～ＬＯＳＪにおいてレーザレーダセンサ１１からの複数の距離範囲Δ（１），Δ（２），…，Δ（Ｎｂｉｎ）にそれぞれ存在する風（流体）の視線方向速度の計測値Ｖｏｂｓ (ｋ)（１，ｊ），…，Ｖｏｂｓ (ｋ)（Ｎｂｉｎ，ｊ）を算出する。風況計測部４４は、視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）と積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）との組合せを風計測データＷＭＤとしてバッファメモリ４５に出力する。バッファメモリ４５は、風計測データＷＭＤを一時的に記憶する。
　次に、風況推定部４６では、レンジビン番号ｉと番号ｊとが初期値「１」に設定される（ステップＳＴ２１，ＳＴ２２）。続いて、有効性評価部６１は、バッファメモリ４５から読み出された風計測データＷＭＤに含まれる積算回数Ｎｓｕｍ (ｋ)（ｉ，ｊ）を用いて、視線方向速度の計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の有効性を評価する（ステップＳＴ２３）。その評価結果である評価値ｑ(ｋ)（ｉ，ｊ）は風速算出部６２に出力される。
　次いで、有効性評価部６１は、番号ｊが上限値Ｊに到達していない場合には（ステップＳＴ２４のＮＯ）、番号ｊを１だけインクリメントし（ステップＳＴ２５）、ステップＳＴ２３を実行する。最終的に、番号ｊが上限値Ｊに到達した場合には（ステップＳＴ２４のＹＥＳ）、番号ｊのすべての値について、計測値Ｖｏｂｓ (ｋ)（ｉ，ｊ）の有効性が評価されている。
　次に、風速算出部６２は、上記のとおり、計測値Ｖｏｂｓ (ｋ)（ｉ，１）～Ｖｏｂｓ (ｋ)（ｉ，Ｊ）及び評価値ｑ(ｋ)（ｉ，１）～ｑ(ｋ)（ｉ，Ｊ）を用いて推定風速値Ｖｃｅｎ (ｋ)（ｉ）を算出する（ステップＳＴ２６）。
　次に、レンジビン番号ｉが上限値Ｎｂｉｎに到達していない場合には（ステップＳＴ２７のＮＯ）、有効性評価部６１は、レンジビン番号ｉを１だけインクリメントし（ステップＳＴ２８）、ステップＳＴ２３～ＳＴ２５を実行する。また風速算出部６２は、ステップＳＴ２６を実行する。最終的に、番号ｉが上限値Ｎｂｉｎに到達した場合には（ステップＳＴ２７のＹＥＳ）、レンジビン番号ｉのすべての値について、推定風速値Ｖｃｅｎ (ｋ)（ｉ）が算出されている。
　次に、データ同化部７４は、データ同化処理を実行する（ステップＳＴ３１）。すなわち、データ同化部７４は、上記のとおり、現時刻ｔｋの予測風速場を、風速推定部６０から供給された現時刻ｔｋの推定風速場に同化させて濾波風速場を算出する。図８は、データ同化処理の詳細な手順の一例を示すフローチャートである。
　図８の例では、データ同化部７４は、先ず、レンジビン番号ｉを「－１」に設定する（ステップＳＴ）。次にデータ同化部７４は、サイクル番号ｋが「１」以下の場合には（ステップＳＴ４２のＮＯ）、ステップＳＴ４３～ＳＴ４９を実行し、サイクル番号ｋが「１」よりも大きい場合には（ステップＳＴ４２のＹＥＳ）、ステップＳＴ５１～ＳＴ５６を実行する。
　以下、サイクル番号ｋが「１」以下の場合（ステップＳＴ４２のＮＯ）について説明する。レンジビン番号ｉが「－１」または「０」の場合は（ステップＳＴ４３のＹＥＳ）、データ同化部７４は、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を「０」に設定する（ステップＳＴ４４）。レンジビン番号ｉが「１」以上の場合であって（ステップＳＴ４３のＮＯ）、推定風速値Ｖｃｅｎ (ｋ)（ｉ）が有効であると評価されている場合には（ステップＳＴ４５のＹＥＳ）、データ同化部７４は、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を推定風速値Ｖｃｅｎ (ｋ)（ｉ）に設定する（ステップＳＴ４６）。
　一方、レンジビン番号ｉが「１」以上の場合であって（ステップＳＴ４３のＮＯ）、推定風速値Ｖｃｅｎ (ｋ)（ｉ）が有効ではないと評価されている場合には（ステップＳＴ４５のＮＯ）、データ同化部７４は、有効な推定風速値Ｖｃｅｎ (ｋ)（ｉｅ）を用いて、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を補間する（ステップＳＴ４６）。たとえば、次式（３７）を用いて濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を補間することができるが、これに限定されるものではない。

Figure JPOXMLDOC01-appb-I000037
　式（３７）において、ｉｕｐは、レンジビン番号ｉの昇順に推定風速値Ｖｃｅｎ (ｋ)（ｉ）を探索したときに最初に発見された有効な推定風速値Ｖｃｅｎ (ｋ)（ｉｕｐ）のレンジビン番号であり、ｉｄｎは、レンジビン番号ｉの降順に推定風速値Ｖｃｅｎ (ｋ)（ｉ）を探索したときに最初に発見された有効な推定風速値Ｖｃｅｎ (ｋ)（ｉｄｎ）のレンジビン番号である。ただし、レンジビン番号ｉｕｐが見つからなかった場合には、推定風速値Ｖｃｅｎ (ｋ)（ｉｕｐ＝Ｎｂｉｎ）が「０」に設定されればよい。
　ステップＳＴ４４，ＳＴ４６，Ｓ４７のいずれかが実行された後、データ同化部７４は、レンジビン番号ｉが上限値Ｎｂｉｎに到達していない場合には（ステップＳＴ４８のＮＯ）、レンジビン番号ｉを１だけインクリメントし（ステップＳＴ４９）、ステップＳＴ４３以後の処理を実行する。最終的に、レンジビン番号ｉが上限値Ｎｂｉｎに到達した場合には（ステップＳＴ４８のＹＥＳ）、データ同化処理は完了する。
　次に、サイクル番号ｋが「１」よりも大きい場合（ステップＳＴ４２のＹＥＳ）について説明する。レンジビン番号ｉが「－１」または「０」の場合は（ステップＳＴ５１のＹＥＳ）、データ同化部７４は、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を、現時刻ｔｋの予測風速値Ｖｅｓｔ (ｋ)（ｉ）に設定する（ステップＳＴ５２）。
　レンジビン番号ｉが「１」以上の場合であって（ステップＳＴ５１のＮＯ）、推定風速値Ｖｃｅｎ (ｋ)（ｉ）が有効ではないと評価されている場合には（ステップＳＴ５３のＮＯ）、データ同化部７４は、濾波風速値Ｖｆｉｌ (ｋ)（ｉ）を、現時刻ｔｋの予測風速値Ｖｅｓｔ (ｋ)（ｉ）に設定する（ステップＳＴ５２）。一方、推定風速値Ｖｃｅｎ (ｋ)（ｉ）が有効であると評価されている場合には（ステップＳＴ５３のＹＥＳ）、データ同化部７４は、上記のとおり、現時刻ｔｋの予測風速値と現時刻ｔｋの推定風速値との重み付け和を濾波風速値Ｖｆｉｌ (ｋ)（ｉ）として生成する（ステップＳＴ５５）。
　ステップＳＴ５２，ＳＴ５５のいずれかが実行された後、データ同化部７４は、レンジビン番号ｉが上限値Ｎｂｉｎに到達していない場合には（ステップＳＴ５６のＮＯ）、レンジビン番号ｉを１だけインクリメントし（ステップＳＴ５７）、ステップＳＴ５１以後の処理を実行する。最終的に、レンジビン番号ｉが上限値Ｎｂｉｎに到達した場合には（ステップＳＴ５６のＹＥＳ）、データ同化処理は完了する。
　図７を参照すると、ステップＳＴ３１のデータ同化処理が完了すると、時間ステップ設定部７７は、ステップＳＴ３１で算出された濾波風速場に基づいて、時間発展計算部７８にて使用される時間ステップｄτを適応的に設定する（ステップＳＴ３２）。次に、時間発展計算部７８は、上記のとおりに、１次元空間に関する時間発展方程式を用いて、濾波風速場からの風速場の時間発展を計算することにより将来時刻ｔｋ＋１，ｔｋ＋Ｃの予測風速場を算出し（ステップＳＴ３３）、将来時刻ｔｋ＋Ｃの予測風速場を予測風速データＰＷＤとして外部に出力する（ステップＳＴ３４）。
　その後、制御部４９は、風況計測を続行する場合には（ステップＳＴ３５のＹＥＳ）、サイクル番号ｋを１だけインクメントし（ステップＳＴ３６）、ステップＳＴ１２に処理を戻す。他方で、風況計測を続行しない場合には（ステップＳＴ３５のＮＯ）、制御部４９は、処理を終了させる。
　以上に説明したように実施の形態１の信号処理器４１は、複数の視線方向ＬＯＳ１～ＬＯＳＪにおける風の視線方向速度の計測値に基づいてレーザレーダセンサ１１の中心軸方向に分布する１次元の推定風速場を算出し、ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、当該１次元の推定風速場から風速場の時間発展を計算する。これにより、高い精度を保ちつつ少ない演算量で予測風速場を算出することが可能となる。
　また、信号処理器４１は低演算負荷で予測風速場を算出することができるので、風力タービンのナセルに搭載可能なコンパクトな構成のレーザレーダ装置１を実現することができる。ナセルに搭載されている制御装置は、レーザレーダ装置１から出力された予測風速データＰＷＤを用いて、ヨー制御及びピッチ制御を行うことが可能である。
　上記した信号処理器４１の機能の全部または一部は、たとえば、ＤＳＰ（Ｄｉｇｉｔａｌ　Ｓｉｇｎａｌ　Ｐｒｏｃｅｓｓｏｒ），ＡＳＩＣ（Ａｐｐｌｉｃａｔｉｏｎ　Ｓｐｅｃｉｆｉｃ　Ｉｎｔｅｇｒａｔｅｄ　Ｃｉｒｃｕｉｔ）またはＦＰＧＡ（Ｆｉｅｌｄ－Ｐｒｏｇｒａｍｍａｂｌｅ　Ｇａｔｅ　Ａｒｒａｙ）などの半導体集積回路を有する単数または複数のプロセッサにより実現可能である。あるいは、信号処理器４１の機能の全部または一部は、ソフトウェアまたはファームウェアのプログラムコードを実行する、ＣＰＵ（Ｃｅｎｔｒａｌ　Ｐｒｏｃｅｓｓｉｎｇ　Ｕｎｉｔ）またはＧＰＵ（Ｇｒａｐｈｉｃｓ　Ｐｒｏｃｅｓｓｉｎｇ　Ｕｎｉｔ）などの演算装置を含む単数または複数のプロセッサで実現されてもよい。あるいは、ＤＳＰ，ＡＳＩＣまたはＦＰＧＡなどの半導体集積回路と、ＣＰＵまたはＧＰＵなどの演算装置との組み合わせを含む単数または複数のプロセッサによって信号処理器４１の機能の全部または一部を実現することも可能である。
　図９は、実施の形態１の信号処理器４１のハードウェア構成例である信号処理回路９０の概略構成を示す機能ブロック図である。図９に示される信号処理回路９０は、プロセッサ９１、入出力インタフェース部９４、メモリ９２、記憶装置９３及び信号路９５を備えている。信号路９５は、プロセッサ９１、入出力インタフェース部９４、メモリ９２及び記憶装置９３を相互に接続するためのバスである。入出力インタフェース部９４は、レーザレーダセンサ１１から入力されたアナログ受信信号Ｒｊ（ｈ，ｔ）をディジタル信号に変換し、当該ディジタル信号をプロセッサ９１に転送する機能を有するとともに、制御信号を外部に出力する機能を有している。
　メモリ９２は、プロセッサ９１がディジタル信号処理を実行する際に使用されるワークメモリと、当該ディジタル信号処理で使用されるデータが展開される一時記憶メモリとを含む。たとえば、メモリ９２は、フラッシュメモリ及びＳＤＲＡＭ（Ｓｙｎｃｈｒｏｎｏｕｓ　Ｄｙｎａｍｉｃ　Ｒａｎｄｏｍ　Ａｃｃｅｓｓ　Ｍｅｍｏｒｙ）などの半導体メモリで構成されればよい。また、記憶装置９３は、プロセッサ９１がＣＰＵまたはＧＰＵなどの演算装置を含む場合には、当該演算装置で実行されるべきソフトウェアまたはファームウェアのプログラムコードを格納する記憶領域として利用可能である。たとえば、記憶装置９３は、フラッシュメモリまたはＲＯＭ（Ｒｅａｄ　Ｏｎｌｙ　Ｍｅｍｏｒｙ）などの不揮発性の半導体メモリで構成されればよい。
　なお、図９の例では、プロセッサ９１の個数は１つであるが、これに限定されるものではない。互いに連携して動作する複数個のプロセッサを用いて信号処理器４１のハードウェア構成が実現されてもよい。
実施の形態２．
　次に、本発明に係る実施の形態２について説明する。上記した実施の形態１では、あらかじめ決められた重み係数α，１－αを用いて濾波風速値Ｖｆｉｌ (ｋ)（ｉ）が算出された。これに対し、以下に説明する実施の形態２，３では、重み係数α，１－αは適応的に決定される。
　図１０は、本発明に係る実施の形態２のレーザレーダ装置における信号処理器４２の概略構成を示すブロック図である。本実施の形態のレーザレーダ装置は、実施の形態１と同様にレーザレーダセンサ１１を備えている。本実施の形態のレーザレーダ装置は、実施の形態１の場合と同様に風力タービンに搭載可能である。
　本実施の形態の信号処理器４２の構成は、図１の風況推定部４６に代えて風況推定部４７を有する点を除いて、実施の形態１の信号処理器４１の構成と同じである。図１１は、本実施の形態における風況推定部４７の概略構成を示すブロック図である。
　図１１に示されるように風況推定部４７は、風速推定部６０及び風速予測部７２を有している。風速予測部７２は、データ同化部７５、時間ステップ設定部７７、時間発展計算部７８及び遅延素子７９を有する。本実施の形態の風況推定部４７の構成は、図４のデータ同化部７４に代えてデータ同化部７５を有する点を除いて、実施の形態１の風況推定部４６の構成と同じである。
　データ同化部７５は、現時刻ｔｋの推定風速場に基づいて、現時刻ｔｋの予測風速場に割り当てるべき重み係数αと現時刻ｔｋの推定風速場に割り当てるべき重み係数１－αとを適応的に設定することを特徴とする。たとえば、データ同化部７５は、現時刻ｔｋの推定風速場を構成する推定風速値Ｖｃｅｎ (ｋ)（ｉ）が大きいほど、重み係数１－αの値を大きくするとともに重み係数αの値を小さくし、現時刻ｔｋの推定風速場を構成する推定風速値Ｖｃｅｎ (ｋ)（ｉ）が小さいほど、重み係数１－αの値を小さくするとともに重み係数αの値を大きくすることができる。これにより、風況に応じて、精度の高い濾波風速場を算出することが可能となる。
　図１２は、実施の形態２に係るデータ同化処理の手順の一例を概略的に示すフローチャートである。図１２のフローチャートは、ステップＳＴ５３とステップＳＴ５５との間にステップＳＴ５４Ａを有する点を除いて、図８のフローチャートと同じである。ステップＳＴ５４Ａでは、データ同化部７５は、現時刻ｔｋの推定風速場に基づいて重み係数αを適応的に設定する。
実施の形態３．
　次に、本発明に係る実施の形態３について説明する。図１３は、本発明に係る実施の形態３のレーザレーダ装置における信号処理器４３の概略構成を示すブロック図である。本実施の形態のレーザレーダ装置は、実施の形態１と同様にレーザレーダセンサ１１を備えている。本実施の形態のレーザレーダ装置は、実施の形態１の場合と同様に風力タービンに搭載可能である。
　本実施の形態の信号処理器４３の構成は、図１の風況推定部４６に代えて風況推定部４８を有する点を除いて、実施の形態１の信号処理器４１の構成と同じである。図１４は、本実施の形態における風況推定部４８の概略構成を示すブロック図である。
　図１４に示されるように風況推定部４８は、風速推定部６０及び風速予測部７３を有している。風速予測部７３は、データ同化部７６、時間ステップ設定部７７、時間発展計算部７８及び遅延素子７９を有する。本実施の形態の風況推定部４８の構成は、図４のデータ同化部７４に代えてデータ同化部７６を有する点を除いて、実施の形態１の風況推定部４６の構成と同じである。
　データ同化部７６は、図１４に示されるように風速推定部６０から、有効な計測値の個数Ｎｑ (ｋ)（ｉ）の供給を受けている。データ同化部７６は、有効な計測値の個数Ｎｑ (ｋ)（ｉ）に基づいて、現時刻ｔｋの予測風速場に割り当てるべき重み係数αと現時刻ｔｋの推定風速場に割り当てるべき重み係数１－αとを適応的に設定することを特徴とする。たとえば、データ同化部７６は、有効な計測値の個数Ｎｑ (ｋ)（ｉ）が多いほど、重み係数１－αの値を大きくするとともに重み係数αの値を小さくし、有効な計測値の個数Ｎｑ (ｋ)（ｉ）が少ないほど、重み係数１－αの値を小さくするとともに重み係数αの値を大きくすることができる。これにより、風況に応じて、精度の高い濾波風速場を算出することが可能となる。
　図１５は、実施の形態３に係るデータ同化処理の手順の一例を概略的に示すフローチャートである。図１５のフローチャートは、ステップＳＴ５３とステップＳＴ５５との間にステップＳＴ５４Ｂを有する点を除いて、図８のフローチャートと同じである。ステップＳＴ５４Ｂでは、データ同化部７６は、有効な計測値の個数Ｎｑ (ｋ)（ｉ）に基づいて重み係数αを適応的に設定する。
実施の形態４．
　次に、本発明に係る実施の形態４について説明する。上記した実施の形態１～３のレーザレーダ装置はいずれも図１に示したレーザレーダセンサ１１を備えているが、これに限定されるものではない。図１６は、本発明に係る実施の形態４のレーザレーダ装置４の概略構成を示すブロック図である。
　図１６に示されるようにレーザレーダ装置４は、レーザレーダセンサ１１とは異なる構成を有するレーザレーダセンサ１２と、信号処理器４１とを備えている。本実施の形態の信号処理器４１の構成は、実施の形態１の信号処理器４１の構成と同じである。本実施の形態の信号処理器４１に代えて、実施の形態２の信号処理器４２または実施の形態３の信号処理器４３のいずれかが使用されてもよい。
　図１６に示されるレーザレーダセンサ１２は、波長可変レーザ光源（ｗａｖｅｌｅｎｇｔｈ－ｔｕｎａｂｌｅ　ｌａｓｅｒ　ｌｉｇｈｔ　ｓｏｕｒｃｅ）２２、光分配器２３、光変調器２４、光増幅器２５、光サーキュレータ２６、光合波器２７、光検出器２８及び光送受信部３２を備えて構成されている。レーザレーダセンサ１２の構成は、図１のレーザ光源２１及び光送受信部３１に代えて、波長可変レーザ光源２２及び光送受信部３２を有する点を除いて、実施の形態１のレーザレーダセンサ１１の構成と同じである。
　本実施の形態の波長可変レーザ光源２２は、信号処理器４１から供給された制御信号に応じて、出力レーザ光の波長λ１，λ２，…，λＪを切替可能な機能を有している。
　光送受信部３２は、光マルチプレクサなどの波長分離素子３４と、光アンテナ３５１～３５Ｊとを備えて構成されている。波長分離素子３４は、入力された増幅光パルスの波長に応じて、光サーキュレータ２６の光入出力端と光アンテナ３５１～３５Ｊのいずれかとの間を光学的に接続する。波長分離素子３４は、ｊ番目の波長λｊを有する増幅光パルスが入力されたときは、当該増幅光パルスをｊ番目の光アンテナ３５ｊにのみ出力する。このとき、当該増幅光パルスは、ｊ番目の光アンテナ３５ｊから視線方向ＬＯＳｊに送出される。観測空間にて当該増幅光パルスが後方散乱することにより散乱光信号が発生すると、ｊ番目の光アンテナ３５ｊは当該散乱光信号を受信光信号として受信する。当該受信光信号は、波長分離素子３４を介して光サーキュレータ２６に伝搬する。
　実施の形態１の光スイッチ３３の切替動作回数には寿命があり、光スイッチ３３よりも波長分離素子３４の寿命の方が長い。レーザ光源２１及び光スイッチ３３に代えて、波長可変レーザ光源２２及び波長分離素子３４を使用することにより、レーザレーダ装置４の長寿命化を実現することができる。
　以上、図面を参照して本発明に係る実施の形態１～４及びこれらの変形例について述べたが、実施の形態１～４及びこれらの変形例は本発明の例示であり、実施の形態１～４以外の様々な実施の形態がありうる。本発明の範囲内において、上記実施の形態１～４の自由な組み合わせ、各実施の形態の任意の構成要素の変形、または各実施の形態の任意の構成要素の省略が可能である。
　たとえば、実施の形態１の場合と同様に、実施の形態２～４の各々の信号処理器の機能の全部または一部は、ＤＳＰ，ＡＳＩＣまたはＰＬＤなどの半導体集積回路を有する単数または複数のプロセッサにより実現可能である。あるいは、信号処理器の機能の全部または一部は、ソフトウェアまたはファームウェアのプログラムコードを実行する、ＣＰＵまたはＧＰＵなどの演算装置を含む単数または複数のプロセッサにより実現されてもよい。あるいは、ＤＳＰ，ＡＳＩＣまたはＰＬＤなどの半導体集積回路と、ＣＰＵまたはＧＰＵなどの演算装置との組み合わせを含む単数または複数のプロセッサによって信号処理器の機能の全部または一部を実現することも可能である。図９に示した信号処理回路９０によって信号処理器のハードウェア構成が実現されてもよい。
　本発明に係る信号処理器，レーザレーダ装置及び風力タービンは、風力タービンを有する風力発電システムまたはウィンドファームが配置された領域における風況観測に利用されることに適している。また本発明に係る信号処理器，レーザレーダ装置及び風力タービンは、風力発電システムまたはウィンドファームの建設予定地における風況観測に利用されることにも適している。
　１，４　レーザレーダ装置、１１，１２　レーザレーダセンサ、２１　レーザ光源、２２　波長可変レーザ光源２２、２３　光分配器、２４　光変調器、２５　光増幅器、２６　光サーキュレータ、２７　光合波器、２８　光検出器、３１，３２　光送受信部、３３　光スイッチ、３４　波長分離素子、３５１～３５Ｊ　光アンテナ、４１，４２，４３　信号処理器、４４　風況計測部、４５　バッファメモリ、４６，４７，４８　風況推定部、４９　制御部、５３　品質評価部、５６　スペクトル積算部、５７　視線速度算出部、６０　風速推定部、６１　有効性評価部、６２　風速算出部、７１，７２，７３　風速予測部、７４，７５，７６　データ同化部、７７　時間ステップ設定部、７８　時間発展計算部、７９　遅延素子、９０　信号処理回路、９１　プロセッサ、９２　メモリ、９３　記憶装置、９４　入出力インタフェース部、９５　信号路、１００　風力タービン、１０１　ナセル、１０２　タワー、１０３　ハブ、１０４Ａ，１０４Ｂ，１０４Ｃ　ブレード、ＬＯＳ１～ＬＯＳＪ　視線方向。

Claims (20)

  1. 　レーザ光源を用いて生成された複数の変調光信号を観測空間の複数方位にそれぞれ送出し、前記複数の変調光信号が後方散乱することにより生じた複数の散乱光信号を複数の受信光信号として受信し、前記複数の受信光信号を検波して前記複数方位にそれぞれ対応する複数の受信信号を生成するレーザレーダセンサと連携して動作する信号処理器であって、
    　前記複数の受信信号を周波数解析することにより、前記複数方位において前記レーザレーダセンサからの複数の距離範囲にそれぞれ存在する風の視線方向速度の計測値を算出する風況計測部と、
    　前記視線方向速度の計測値に基づいて、前記レーザレーダセンサの中心軸方向に分布する１次元風速場を現時刻の推定風速場として算出する風速推定部と、
    　ナビエ・ストークス方程式から導出された１次元空間に関する時間発展方程式を用いて、前記推定風速場からの風速場の時間発展を計算することにより将来時刻の予測風速場を算出する風速予測部と
    を備えることを特徴とする信号処理器。
  2. 　請求項１に記載の信号処理器であって、前記時間発展方程式は、前記風が非圧縮流体であり、前記風の粘性が考慮されず、かつ外力が存在しないとの条件に基づいて前記ナビエ・ストークス方程式から導出された方程式であることを特徴とする信号処理器。
  3. 　請求項２に記載の信号処理器であって、前記時間発展方程式は、前記中心軸方向とは垂直な径方向における距離に関する一階偏微分が考慮されないとの条件に基づいて前記ナビエ・ストークス方程式から導出された方程式であることを特徴とする信号処理器。
  4. 　請求項３に記載の信号処理器であって、前記時間発展方程式は、前記距離に関する圧力の二階偏微分が当該圧力に関する線形関数であるとの条件に基づいて前記ナビエ・ストークス方程式から導出された方程式であることを特徴とする信号処理器。
  5. 　請求項２から請求項４のうちのいずれか１項に記載の信号処理器であって、前記時間発展方程式は、風上差分スキームに従って前記１次元空間において離散化された式を含むことを特徴とする信号処理器。
  6. 　請求項１から請求項５のうちのいずれか１項に記載の信号処理器であって、
    　前記風速予測部は、
    　当該風速予測部により以前に算出された現時刻の予測風速場を、前記風速推定部により算出された現時刻の推定風速場に同化させて濾波風速場を生成するデータ同化部と、
    　前記濾波風速場からの風速場の時間発展を計算することにより前記将来時刻の予測風速場を算出する時間発展計算部と
    を含むことを特徴とする信号処理器。
  7. 　請求項６に記載の信号処理器であって、
    　前記データ同化部は、前記現時刻の予測風速場と前記現時刻の推定風速場との重み付け和を前記濾波風速場として生成することを特徴とする信号処理器。
  8. 　請求項６または請求項７に記載の信号処理器であって、
    　前記風速予測部は、前記濾波風速場に基づいて時間ステップを適応的に設定する時間ステップ設定部をさらに含み、
    　前記時間発展計算部は、前記時間ステップにて前記濾波風速場からの風速場の時間発展を計算することにより前記将来時刻の予測風速場を算出する、
    ことを特徴とする信号処理器。
  9. 　請求項８に記載の信号処理器であって、前記時間ステップ設定部は、前記濾波風速場を構成する複数の風速値のうちの最大風速値に基づいて前記時間ステップを適応的に設定することを特徴とする信号処理器。
  10. 　請求項１から請求項９のうちのいずれか１項に記載の信号処理器であって、前記風速推定部は、前記視線方向速度の計測値を用いた非線形最小自乗法により前記推定風速場を算出することを特徴とする信号処理器。
  11. 　請求項１０に記載の信号処理器であって、前記風速推定部は、前記複数の距離範囲にそれぞれ存在する当該風の速度ベクトルの方向が前記中心軸方向と平行であるとの仮定の下で、前記非線形最小自乗法により前記推定風速場を算出することを特徴とする信号処理器。
  12. 　請求項１から請求項１１のうちのいずれか１項に記載の信号処理器であって、
    　前記風速推定部は、
    　前記視線方向速度の計測値それぞれの有効性を評価して評価結果を出力する有効性評価部と、
    　前記評価結果に基づき、前記視線方向速度の計測値のうち前記複数の距離範囲の各距離範囲に対応する有効な計測値の個数が閾値以上であるか否かを判定する風速算出部と
    を含み、
    　前記風速算出部は、前記有効な計測値の個数が閾値以上であると判定したときは、前記有効な計測値に基づいて推定風速値を算出し、前記有効な計測値の個数が閾値未満であると判定したときは、推定風速値を算出せず、
    　前記推定風速場は、前記複数の距離範囲について前記風速算出部により算出された単数または複数の推定風速値からなる、
    ことを特徴とする信号処理器。
  13. 　請求項７に記載の信号処理器であって、前記データ同化部は、前記現時刻の推定風速場に基づいて、前記現時刻の予測風速場に割り当てるべき第１の重み係数と前記現時刻の推定風速場に割り当てるべき第２の重み係数とを適応的に設定することを特徴とする信号処理器。
  14. 　請求項１３に記載の信号処理器であって、前記データ同化部は、前記現時刻の推定風速場を構成する推定風速値が大きいほど、前記第２の重み係数の値を大きくするとともに前記第１の重み係数の値を小さくし、前記現時刻の推定風速場を構成する推定風速値が小さいほど、前記第２の重み係数の値を小さくするとともに前記第１の重み係数の値を大きくする、ことを特徴とする信号処理器。
  15. 　請求項７に記載の信号処理器であって、
    　前記風速推定部は、
    　前記視線方向速度の計測値それぞれの有効性を評価して評価結果を出力する有効性評価部と、
    　前記評価結果に基づき、前記視線方向速度の計測値のうち前記複数の距離範囲の各距離範囲に対応する有効な計測値の個数が閾値以上であるか否かを判定する風速算出部と
    を含み、
    　前記風速算出部は、前記有効な計測値の個数が閾値以上であると判定したときは、前記有効な計測値に基づいて推定風速値を算出し、前記有効な計測値の個数が閾値未満であると判定したときは、推定風速値を算出せず、
    　前記推定風速場は、前記複数の距離範囲に対して前記風速算出部により算出された単数または複数の推定風速値からなり、
    　前記データ同化部は、前記単数または複数の風速値の個数に基づいて、前記現時刻の予測風速場に割り当てるべき第１の重み係数と前記現時刻の推定風速場に割り当てるべき第２の重み係数とを適応的に設定する、ことを特徴とする信号処理器。
  16. 　請求項１５に記載の信号処理器であって、前記データ同化部は、前記単数または複数の風速値の個数が多いほど、前記第２の重み係数の値を大きくするとともに前記第１の重み係数の値を小さくし、前記単数または複数の風速値の個数が少ないほど、前記第２の重み係数の値を小さくするとともに前記第１の重み係数の値を大きくする、ことを特徴とする信号処理器。
  17. 　請求項１２または請求項１５に記載の信号処理器であって、
    　前記風況計測部は、
    　前記複数の受信信号の各受信信号を、前記複数の距離範囲にそれぞれ対応する複数の時間領域信号に分割する時間ゲーティング部と、
    　前記複数の時間領域信号に対して時間領域から周波数領域への領域変換処理を実行することにより複数の周波数領域信号を生成する領域変換部と、
    　前記複数の周波数領域信号それぞれの信号品質を評価する品質評価部と
    を含み、
    　前記有効性評価部は、前記品質評価部により評価された信号品質に基づいて前記視線方向速度の計測値それぞれの有効性を評価する、ことを特徴とする信号処理器。
  18. 　請求項１から請求項１７のうちのいずれか１項に記載の信号処理器と、
    　前記レーザレーダセンサと
    を備えることを特徴とするレーザレーダ装置。
  19. 　請求項１８に記載のレーザレーダ装置と、
    　前記レーザレーダ装置を搭載するナセルと、
    　前記ナセルに取り付けられたハブと、
    　前記ハブを介して前記ナセルと回転自在に連結された複数枚のブレードからなるロータと
    を備えることを特徴とする風力タービン。
  20. 　請求項１９に記載の風力タービンであって、前記レーザレーダ装置は、前記ロータの回転面を介して前記観測空間に前記複数の変調光信号を送出するように前記ナセル上に配置されていることを特徴とする風力タービン。
PCT/JP2019/042018 2019-10-25 2019-10-25 信号処理器、レーザレーダ装置及び風力タービン WO2021079513A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
PCT/JP2019/042018 WO2021079513A1 (ja) 2019-10-25 2019-10-25 信号処理器、レーザレーダ装置及び風力タービン
JP2021553876A JP6995256B2 (ja) 2019-10-25 2019-10-25 信号処理器、レーザレーダ装置及び風力タービン

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2019/042018 WO2021079513A1 (ja) 2019-10-25 2019-10-25 信号処理器、レーザレーダ装置及び風力タービン

Publications (1)

Publication Number Publication Date
WO2021079513A1 true WO2021079513A1 (ja) 2021-04-29

Family

ID=75619758

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2019/042018 WO2021079513A1 (ja) 2019-10-25 2019-10-25 信号処理器、レーザレーダ装置及び風力タービン

Country Status (2)

Country Link
JP (1) JP6995256B2 (ja)
WO (1) WO2021079513A1 (ja)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114167522A (zh) * 2021-12-17 2022-03-11 深圳市云端高科信息科技有限公司 一种用于智慧城市高层建筑的风场监测与校正系统
JP7069433B1 (ja) * 2021-06-21 2022-05-17 三菱電機株式会社 風速予測装置、風速予測方法及びレーダ装置
CN116540266A (zh) * 2023-07-06 2023-08-04 南京牧镭激光科技股份有限公司 基于相干激光雷达识别极端风况的方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2004285858A (ja) * 2003-03-19 2004-10-14 Mitsubishi Electric Corp 風力発電システムおよび風力発電機の制御方法
JP2009052961A (ja) * 2007-08-24 2009-03-12 Mitsubishi Electric Corp 風計測装置
US20110216307A1 (en) * 2010-02-05 2011-09-08 Catch the Wind, Inc. High Density Wind Velocity Data Collection for Wind Turbine
EP2581761A1 (en) * 2011-10-14 2013-04-17 Vestas Wind Systems A/S Estimation of Wind Properties Using a Light Detection and Ranging Device

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2004285858A (ja) * 2003-03-19 2004-10-14 Mitsubishi Electric Corp 風力発電システムおよび風力発電機の制御方法
JP2009052961A (ja) * 2007-08-24 2009-03-12 Mitsubishi Electric Corp 風計測装置
US20110216307A1 (en) * 2010-02-05 2011-09-08 Catch the Wind, Inc. High Density Wind Velocity Data Collection for Wind Turbine
EP2581761A1 (en) * 2011-10-14 2013-04-17 Vestas Wind Systems A/S Estimation of Wind Properties Using a Light Detection and Ranging Device

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
RAACH, STEFFEN ET AL.: "Three Dimensional Dynamic Model BasedWind Field Reconstruction from Lidar Data", JOURNAL OF PHYSICS: CONFERENCE SERIES 524 (TORQUE 2014), vol. 524, no. 1, 18 June 2014 (2014-06-18) - 20 June 2014 (2014-06-20), pages 012005, XP020265679, DOI: 10.1088/1742- 6596/524/1/012005 *
RETHORE, P. E. ET AL.: "A discrete force allocation algorithm for modelling wind turbines in computational fluid dynamics", WIND ENERGY, vol. 15, no. 7, October 2012 (2012-10-01), pages 915 - 926, XP055819722, DOI: 10.1002/we.525 *
TOWERS, P. ET AL.: "Real-time wind field reconstruction from LiDAR measurements using a dynamic wind model and state estimation", WIND ENERGY, vol. 19, no. 1, 2016, pages 133 - 150, XP055451329, DOI: 10.1002/we.1824 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP7069433B1 (ja) * 2021-06-21 2022-05-17 三菱電機株式会社 風速予測装置、風速予測方法及びレーダ装置
WO2022269655A1 (ja) * 2021-06-21 2022-12-29 三菱電機株式会社 風速予測装置、風速予測方法及びレーダ装置
CN114167522A (zh) * 2021-12-17 2022-03-11 深圳市云端高科信息科技有限公司 一种用于智慧城市高层建筑的风场监测与校正系统
CN116540266A (zh) * 2023-07-06 2023-08-04 南京牧镭激光科技股份有限公司 基于相干激光雷达识别极端风况的方法
CN116540266B (zh) * 2023-07-06 2023-10-27 南京牧镭激光科技股份有限公司 基于相干激光雷达识别极端风况的方法

Also Published As

Publication number Publication date
JP6995256B2 (ja) 2022-01-14
JPWO2021079513A1 (ja) 2021-04-29

Similar Documents

Publication Publication Date Title
JP6995256B2 (ja) 信号処理器、レーザレーダ装置及び風力タービン
KR102208288B1 (ko) Lidar 센서에 기초한 풍속 추정을 사용한 풍력 터빈 제어 및 모니터링 방법
Peña et al. Turbulence characterization from a forward-looking nacelle lidar
Banner et al. Wave breaking onset and strength for two-dimensional deep-water wave groups
Raach et al. Lidar-based wake tracking for closed-loop wind farm control
Wilson Performance bounds for acoustic direction-of-arrival arrays operating in atmospheric turbulence
CN104101731A (zh) 用于远程风感测的方法、装置和系统
ES2882299T3 (es) Procedimiento de determinación de un factor de inducción para un aerogenerador de un sensor de detección a distancia por láser
CN104792364A (zh) 基于激光多普勒的桥梁动态参数提取系统及方法
Leuenberger et al. Revisiting the latent heat nudging scheme for the rainfall assimilation of a simulated convective storm
ES2970376T3 (es) Método de determinación de la velocidad media del viento mediante un sensor de teledetección por láser
Frehlich et al. Measurements of wind and turbulence profiles with scanning Doppler lidar for wind energy applications
Gräfe et al. Wind field reconstruction using nacelle based lidar measurements for floating wind turbines
CN114239109A (zh) 基于节段模型测振试验的大跨度桥梁抖振响应直接预测方法及系统、存储介质
Cheong et al. Cyclostationary spectral analysis for the measurement and prediction of wind turbine swishing noise
Gräfe et al. Quantification and correction of motion influence for nacelle-based lidar systems on floating wind turbines
Zhang et al. Spatial estimation of unidirectional wave evolution based on ensemble data assimilation
CN115469378B (zh) 一种湿位涡测量方法、装置、设备及介质
REN et al. Effect of peak power and pulse width on coherent Doppler wind lidar’s SNR
JP7459408B2 (ja) レーザレーダ装置
CN110290454A (zh) 基于光学法的传声器高温校准系统
CN114839646A (zh) 一种基于激光多普勒联合闪烁的测风装置以及测风方法
Castellani et al. On the way to harness high-altitude wind power: Defining the operational asset for an airship wind generator
Zhang et al. Novel rotor effective wind speed estimation method for light detection and ranging application
CN113063961A (zh) 一种超声波传感阵列测风装置及其方法

Legal Events

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

Ref document number: 19949654

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2021553876

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 19949654

Country of ref document: EP

Kind code of ref document: A1