Bonjour,

J'ai hésité à placer ce post dans la section 'Calcul Scientique', mais je pense que le problème est plus général.

L'intitulé n'est pas très explicite... Je cherche à générer une 'time serie' (ou encore un signal temporel) polychromatique (plusieurs fréquences) comme ceci:

Code : Sélectionner tout - Visualiser dans une fenêtre à part
F(t) = Somme(Ai*cos(wi*t+Phii))
J'ai un certain nombre de "composantes" à mon signal (les indices i) sur une certaine plage de temps. Le résutat est sous forme de liste (un élément par pas de temps t).

Mon problème vient du temps d'éxécution. Beaucoup trop long je trouve (une dizaine de secondes quand même), je pense donc que ma méthode n'est pas optimale.

Voici une version épurée de mon code:
Code : Sélectionner tout - Visualiser dans une fenêtre à part
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
# -*- coding:Utf-8 -*-
 
 
#-----------------------------------------------------------------------
#			        MODULE IMPORTATION
#-----------------------------------------------------------------------
 
import numpy as np
from math import cos
import random
import time
 
#-----------------------------------------------------------------------
#			        SPECTRUM DEFINITIONS
#-----------------------------------------------------------------------
 
 
def HarrisDNV(meanSpeed, Wmin, Wmax, N):
    """ Harris - DNV Spectrum """
 
    # Frequencies in rad/s
    W = linspace(Wmin,Wmax,N)
 
    # Surface drag coefficient  	
    C = 0.0015 
 
    # Spectrum					
    Spectrum = [(7200/(2*np.pi))*C*meanSpeed*(2+(286.0*W[i]/meanSpeed)**2)**(-5.0/6.0) for i in range(len(W))]
 
    return W,Spectrum
 
 
 
#-----------------------------------------------------------------------
#			       RANDOM TIME SERIE GENERATOR
#-----------------------------------------------------------------------
 
def WindGen(Npoints, SampleTime, MaxFreq, Nlags, meanSpeed, Spec=2, Seed=None):
    """ Build a time serie of the wind velocity 
    
    This function generate a random time serie of wind speed based on 
    the following wind spectrum:
        1- Ochi-Shin
        2- Harris-DNV (default)
        3- Wills or Modified Harris
        4- API
        5- NPD
    
    --------------------------------------------------------------------
    INPUTS:
    --------------------------------------------------------------------
        Npoints     : Number of points in the time serie (integer)
        SampleTime  : Sample time (sec) (integer)
        MaxFreq     : Max frequency for spectrum (float)
        Nlags       : Number of random phases (integer)
        meanSpeed   : Mean wind velocity (float)
        Spec        : Spectrum type (integer)
        
    --------------------------------------------------------------------
    OUTPUTS
    --------------------------------------------------------------------
        liste containing:
            - time
            - wind velocities
    """
 
    #
    # GENERATING WIND SPECTRUM
    #
    if Spec == 1:
        print 'OchiShin'
        W,S = OchiShin(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
    elif Spec == 2:
        print 'HarrisDNV'
        W,S = HarrisDNV(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
    elif Spec == 3:
        print 'Wills'
        W,S = Wills(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
    elif Spec == 4:
        print 'API'
        W,S = API(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
    elif Spec == 5:
        print 'NPD'
        W,S = NPD(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
    else:
        print 'Unvalide choice, default spectrum have been selected (Harris)'
        print 'The following spectrum are available:'
        print '\t' + '1- Ochi-Shin'
        print '\t' + '2- Harris-DNV (default)'
        print '\t' + '3- Wills or Modified Harris'
        print '\t' + '4- API'
        print '\t' + '5- NPD'
        print ''
        W,S = HarrisDNV(meanSpeed, Wmin=0., Wmax=MaxFreq, N=Nlags)
 
    #           
    # GENERATING TIME SERIES
    #
    if Seed: random.seed(Seed)
 
    T = [t for t in range(Npoints)]
    A = [(2*(S[i+1]+S[i])/2*(W[i+1]-W[i]))**0.5 for i in range(Nlags-1)]
    A.append(0.0)
    P = [random.uniform(0.0,2*np.pi) for i in range(Nlags)]
 
    t0 = time.clock()
    WindVelocities = [sum([A[i]*cos(W[i]*t + P[i]) for i in range(Nlags)]) + meanSpeed for t in range(Npoints)]
    print 'Time serie generated in ' + str(round(time.clock()-t0,4)) + 's'
 
    return T, WindVelocities
 
 
 
if __name__=="__main__":
 
    import numpy as np
    from pylab import *	
    from math import cos
    import random
    import time
 
    t, WindVelocities = WindGen(Npoints = 3600, 
                                   SampleTime = 1, 
                                   MaxFreq = 2*np.pi/4, 
                                   Nlags = 1000, 
                                   meanSpeed = 1.0, 
                                   Spec=2, 
                                   Seed=145468)
 
 
    plot(t[100:300],WindVelocities[100:300],'blue')
    legend(('Random Wind velocities'),'best')
    ylabel('Wind Velocity in m/s')
    xlabel('Time in s')
    tight_layout()
    show()

La création de la liste 'WindVelocities' prend chez moi plus de 10s... Est-ce que la liste compréhension n'est pas adapté, ou est-ce la fonction 'sum' qui prend autant de temps?



Ju