2011-05-24 97 views
9

我被要求使用rainflow counting algorithm作爲我正在進行的項目。有沒有人有任何我可以用來完成這項任務的庫的經驗。我已經使用Python和Numpy計算了一些其他統計數據,因此使用/擴展這些工具的任何建議都是理想的。雨流計數算法

謝謝。

+2

+1我喜歡被通過StackOverflow的接觸到新的概念。謝謝 – I82Much 2011-05-24 10:14:37

+0

我很想有一個很好的雨點流量計算算法。 – 2011-05-24 19:21:04

回答

4

快速搜索打開了以下內容:

WAFO:工具箱進行統計分析和隨機波和隨機負荷的模擬。這包括雨流計數算法

http://code.google.com/p/pywafo/

+0

這看起來非常好。 Python仍然需要Matlab來使用它嗎? Matlab版本看起來很好記錄,但我只能找到Python版本的SVN代碼,並且我不太清楚應該如何使用它。 – dwxw 2011-05-24 15:16:59

+0

@dww我沒有看到python版本的工具包的matlab依賴關係 – JoshAdel 2011-05-24 17:17:46

7

我不得不安裝WAFO(防火牆有心計)問題的實現,所以我用C移植版本Matlab的雨流計數算法到Python。它運行Endo算法,也可以選擇使用Goodman校正。只有依賴關係是numpy> = 1.3。簡單但可能有用。

Python的雨流功能(在Gist的功能和可用的演示腳本ZIP):

#!/usr/bin/env python 
""" 
------------------------------------------------------------------------------- 
Rainflow counting function with Goodman correction 
Copyright (C) 2015 Jennifer Rinker 
This program is free software: you can redistribute it and/or modify 
it under the terms of the GNU General Public License as published by 
the Free Software Foundation, either version 3 of the License, or 
(at your option) any later version. 
This program is distributed in the hope that it will be useful, 
but WITHOUT ANY WARRANTY; without even the implied warranty of 
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 
GNU General Public License for more details. 
You should have received a copy of the GNU General Public License 
along with this program. If not, see <http://www.gnu.org/licenses/>. 
Contact: Jennifer Rinker, Duke University 
Email: jennifer.rinker-at-duke.edu 
------------------------------------------------------------------------------- 
USAGE: 
To call the function in a script on array of turning points <array_ext>: 
    import rainflow as rf 
    array_out = rf.rainflow(array_ext) 
To run the demonstration from a Python console: 
    >>> execfile('demo_rainflow.py') 
From the terminal (Windows or UNIX-based): 
    $ python demo_rainflow.py 
------------------------------------------------------------------------------- 
DEPENDENCIES: 
- Numpy >= v1.3 
------------------------------------------------------------------------------- 
NOTES: 
Python code modified from rainflow.c code with mex function for Matlab from 
WISDEM project: https://github.com/WISDEM/AeroelasticSE/tree/master/src/AeroelasticSE/rainflow 

Original c code notes: 
/* RAINFLOW $ Revision: 1.1 $ */ 
/* by Adam Nieslony, 2009  */ 
/* The original code has been modified by Gregory Hayman 2012 as follows: */ 
/* abs() has been replaced everywhere with fabs()      */ 
/* the function now applies the Goodman correction to the damage cycle */ 
/*  load ranges using a user-supplied fixed load mean, or a fixed load */ 
/*  mean of zero.              */ 
/* the user can supply a the value of a partial damage cycle: uc_mult */ 
------------------------------------------------------------------------------- 
""" 
from numpy import fabs as fabs 
import numpy as np 

def rainflow(array_ext, 
      flm=0, l_ult=1e16, uc_mult=0.5): 
    """ Rainflow counting of a signal's turning points with Goodman correction 

     Args: 
      array_ext (numpy.ndarray): array of turning points 

     Keyword Args: 
      flm (float): fixed-load mean [opt, default=0] 
      l_ult (float): ultimate load [opt, default=1e16] 
      uc_mult (float): partial-load scaling [opt, default=0.5] 

     Returns: 
      array_out (numpy.ndarray): (5 x n_cycle) array of rainflow values: 
             1) load range 
             2) range mean 
             3) Goodman-adjusted range 
             4) cycle count 
             5) Goodman-adjusted range with flm = 0 

    """ 

    flmargin = l_ult - fabs(flm)   # fixed load margin 
    tot_num = array_ext.size    # total size of input array 
    array_out = np.zeros((5, tot_num-1)) # initialize output array 

    pr = 0         # index of input array 
    po = 0         # index of output array 
    j = -1         # index of temporary array "a" 
    a = np.empty(array_ext.shape)   # temporary array for algorithm 

    # loop through each turning point stored in input array 
    for i in range(tot_num): 

     j += 1     # increment "a" counter 
     a[j] = array_ext[pr] # put turning point into temporary array 
     pr += 1     # increment input array pointer 

     while ((j >= 2) & (fabs(a[j-1] - a[j-2]) <= \ 
       fabs(a[j] - a[j-1]))): 
      lrange = fabs(a[j-1] - a[j-2]) 

      # partial range 
      if j == 2: 
       mean  = (a[0] + a[1])/2. 
       adj_range = lrange * flmargin/(l_ult - fabs(mean)) 
       adj_zero_mean_range = lrange * l_ult/(l_ult - fabs(mean)) 
       a[0]=a[1] 
       a[1]=a[2] 
       j=1 
       if (lrange > 0): 
        array_out[0,po] = lrange 
        array_out[1,po] = mean 
        array_out[2,po] = adj_range 
        array_out[3,po] = uc_mult 
        array_out[4,po] = adj_zero_mean_range 
        po += 1 

      # full range 
      else: 
       mean  = (a[j-1] + a[j-2])/2. 
       adj_range = lrange * flmargin/(l_ult - fabs(mean)) 
       adj_zero_mean_range = lrange * l_ult/(l_ult - fabs(mean)) 
       a[j-2]=a[j] 
       j=j-2 
       if (lrange > 0): 
        array_out[0,po] = lrange 
        array_out[1,po] = mean 
        array_out[2,po] = adj_range 
        array_out[3,po] = 1.00 
        array_out[4,po] = adj_zero_mean_range 
        po += 1 

    # partial range 
    for i in range(j): 
     lrange = fabs(a[i] - a[i+1]); 
     mean  = (a[i] + a[i+1])/2. 
     adj_range = lrange * flmargin/(l_ult - fabs(mean)) 
     adj_zero_mean_range = lrange * l_ult/(l_ult - fabs(mean)) 
     if (lrange > 0): 
      array_out[0,po] = lrange 
      array_out[1,po] = mean 
      array_out[2,po] = adj_range 
      array_out[3,po] = uc_mult 
      array_out[4,po] = adj_zero_mean_range 
      po += 1 

    # get rid of unused entries 
    array_out = array_out[:,:po] 

    return array_out