Télécharger dspr.eso

Retour à la liste

Numérotation des lignes :

  1. C DSPR SOURCE BP208322 16/11/18 21:16:32 9177
  2. SUBROUTINE DSPR
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8(A-H,O-V)
  5. IMPLICIT COMPLEX*16 (W-Z)
  6. CHARACTER *72 TI
  7. CHARACTER*12 MOTX,MOTY
  8. CHARACTER*4 MCLE(2)
  9. LOGICAL INV
  10. C
  11. C=======================================================================
  12. C = CALCUL DE LA DENSITE SPECTRALE DE PUISSANCE D'UN SIGNAL =
  13. C = =
  14. C = SYNTAXE : =
  15. C = =
  16. C = SPEC = DSPR EXP2 EVO1 FMIN V1 FMAX V2 (COUL) ; =
  17. C = =
  18. C = =
  19. C = EXP2 : EXPOSANT DANS NPOINT=2^EXP2, NPOINT ETANT =
  20. C = LE NOMBRE DE REELS DANS LISTREEL =
  21. C = EVO1 : OBJET DE TYPE EVOLUTIO CONTENANT LE SIGNAL A TRAITER=
  22. C = ( UNE COURBE SEULEMENT ) =
  23. C = FMIN : MOT-CLE =
  24. C = V1 : FREQUENCE MINIE A VISUALISER =
  25. C = FMAX : MOT-CLE =
  26. C = V2 : FREQUENCE MAXIE A VISUALISER =
  27. C = COUL : COULEUR ATTRIBUEE A L'OBJET CREE (FACULTATIF) =
  28. C = =
  29. C = =
  30. C = CREATION : 01/04/87 =
  31. C = PROGRAMMEUR : BEAUFILS =
  32. C = MODIFICATIO : PEG WARNING LE IRET2.... 22/2/90 =
  33. C = MODIFICATIO : PEG FREQUENCES CORRECTES 1/6/90 =
  34. C = MODIFICATIO : PEG PERIODOGRAM CORRECT 1/6/90 =
  35. C=======================================================================
  36. C
  37. -INC CCGEOME
  38. -INC CCOPTIO
  39. -INC SMEVOLL
  40. -INC SMLREEL
  41. C
  42. SEGMENT MTRAV1
  43. IMPLIED XX(NCOU),YY(NCOU)
  44. ENDSEGMENT
  45. C
  46. SEGMENT MTRAV2
  47. IMPLIED W(NEXP)
  48. ENDSEGMENT
  49. C
  50. DATA MOTX/'FREQ HZ'/
  51. DATA MOTY/'DSP UN**2/HZ'/
  52. DATA MCLE(1),MCLE(2)/'FMIN','FMAX'/
  53. C
  54. C LECTURE DE EXP2
  55. C
  56. CALL LIRENT(N2,1,IRET1)
  57. IF(IRET1.EQ.0)GOTO 666
  58. C
  59. C LECTURE DE L'OBJET EVOLUTIO CONTENANT LE SIGNAL
  60. C
  61. CALL LIROBJ('EVOLUTIO',IPSIG,1,IRET2)
  62. C PEG IF(IRET1.EQ.0)GOTO 666
  63. IF(IRET2.EQ.0)GOTO 666
  64. C
  65. C LECTURE DE LA FREQUENCE MINIE
  66. C
  67. CALL LIRMOT(MCLE(1),1,IRET,0)
  68. IF(IRET.EQ.1) THEN
  69. CALL LIRREE(FRMI,0,IRET1)
  70. IF(IRET1.EQ.0) THEN
  71. MOTERR(1:4)=MCLE(1)
  72. CALL ERREUR(166)
  73. GOTO 666
  74. ENDIF
  75. IF(FRMI.LT.(-1.D-20)) FRMI=-FRMI
  76. IF(ABS(FRMI).LT.1.D-20) FRMI=-1.D0
  77. ELSE
  78. FRMI=-1.D0
  79. ENDIF
  80. C
  81. C LECTURE DE LA FREQUENCE MAXIE
  82. C
  83. CALL LIRMOT(MCLE(2),1,IRET,0)
  84. IF(IRET.EQ.1) THEN
  85. CALL LIRREE(FRMA,0,IRET1)
  86. IF(IRET1.EQ.0) THEN
  87. MOTERR(1:4)=MCLE(2)
  88. CALL ERREUR(166)
  89. GOTO 666
  90. ENDIF
  91. IF(FRMA.LT.(-1.D-20)) FRMA=-FRMA
  92. IF(ABS(FRMA).LT.1.D-20) FRMA=-1.D0
  93. ELSE
  94. FRMA=-1.D0
  95. ENDIF
  96. C
  97. C LECTURE DE LA COULEUR
  98. C
  99. CALL LIRMOT(NCOUL,NBCOUL,ICOUL,0)
  100. IF(ICOUL.EQ.0) ICOUL=IDCOUL+1
  101. ICOUL=ICOUL-1
  102. C
  103. IF(IERR.NE.0) GOTO 666
  104. C
  105. MEVOL1=IPSIG
  106. SEGACT MEVOL1
  107. KEVOL1=MEVOL1.IEVOLL(1)
  108. SEGACT KEVOL1
  109. MLREE1=KEVOL1.IPROGX
  110. MLREE2=KEVOL1.IPROGY
  111. SEGACT MLREE1,MLREE2
  112. C
  113. L1=MLREE1.PROG(/1)
  114. DT=MLREE1.PROG(2)-MLREE1.PROG(1)
  115. SEGDES MLREE1
  116. C
  117. NPOINT=2**N2
  118. IF(NPOINT.GT.L1) THEN
  119. IF(IIMPI.EQ.1) WRITE(IOIMP,1000) L1,N2,NPOINT
  120. 1000 FORMAT(1H ,'LE NOMBRE DE POINTS DANS LISTEMPS ',I6, ' EST ',
  121. & 'INFERIEURE @ 2**',I5,
  122. & /' ON PRENDRA UNE LONGUEUR DE LISTEMPS DE ',I6,' MOTS ',
  123. & /' ET ON COMPLETERA PAR DES ZEROS')
  124. ELSE
  125. IF(NPOINT.LT.L1) THEN
  126. IF(IIMPI.EQ.1) WRITE(IOIMP,1010) N2
  127. 1010 FORMAT(1H ,'ON TRONQUE LE SIGNAL A 2**',I5,' MOTS',/)
  128. ELSE
  129. IF(IIMPI.EQ.1) WRITE(IOIMP,1030)N2,NPOINT
  130. 1030 FORMAT(1H ,'LA LONGUEUR DU LISTEMP EST EGALE A 2**',I5,
  131. & ' = ',I6,/)
  132. ENDIF
  133. ENDIF
  134. C
  135. NCOU=NPOINT
  136. SEGINI MTRAV1
  137. C
  138. NEXP=NPOINT/2
  139. SEGINI MTRAV2
  140. C
  141. C CHARGEMENT DES TABLEAUX DE TRAVAIL
  142. C
  143. IND1=L1
  144. IF(NPOINT.LT.L1)IND1 = NPOINT
  145. DO 10 I=1,IND1
  146. XX(I)=MLREE2.PROG(I)
  147. C IF(IIMPI.EQ.1)WRITE(IOIMP,*) ' XX(',I,') = ',XX(I)
  148. 10 CONTINUE
  149. IF(NPOINT.GT.L1) THEN
  150. L2=L1+1
  151. DO 11 I=L2,NPOINT
  152. XX(I)=0.D0
  153. 11 CONTINUE
  154. ENDIF
  155. C
  156. C CALCUL DE LA FFT
  157. C
  158. IF(IIMPI.EQ.1) WRITE(IOIMP,*)' CALCUL DE LA FFT DU SIGNAL '
  159. INV=.FALSE.
  160. CALL WEXP(INV,NPOINT,W(1))
  161. CALL FFTL(XX(1),YY(1),W(1),NPOINT)
  162. C
  163. IF(IIMPI.EQ.1) WRITE(IOIMP,*)' FFT CALCULEE '
  164. SEGDES MLREE2
  165. C
  166. C CREATION ET CALCUL DES LISTES DE LA DSP
  167. C
  168. IF(IIMPI.EQ.1) WRITE(IOIMP,*)' CALCUL DE LA DSP DU SIGNAL '
  169. DUREE=DT*DBLE(NPOINT)
  170. DFREQ=1.D0/DUREE
  171. KHALF=NEXP+1
  172. CPP +
  173. KHALF1=KHALF
  174. CPP +
  175. KMIN=INT(FRMI/DFREQ)+1
  176. KMAX=INT(FRMA/DFREQ)+1
  177. KDEBU=1
  178. IF(KMIN.GT.0) KDEBU=KMIN
  179. IF((KMAX.LT.KHALF).AND.(FRMA.GT.0.)) KHALF=KMAX
  180. JG=KHALF-KDEBU+1
  181. SEGINI MLREE1
  182. SEGINI MLREE2
  183. CPP ! COEF=2.D0*DFREQ*DT*DT
  184. COEF= DFREQ*DT*DT
  185. DO 20 I=KDEBU,KHALF
  186. FREQ=DBLE(I-1)*DFREQ
  187. C IF(IIMPI.EQ.1) WRITE(IOIMP,*)I,FREQ,' XX = ',XX(I)
  188. MLREE1.PROG(I-KDEBU+1)=FREQ
  189. YCOMP=CONJG(XX(I))
  190. PMODU2=XX(I)*YCOMP
  191. cbp: on pourrait aussi ecrire PMODU2=ABS(XX(I))**2
  192. CPP +
  193. IF(I.NE.1.AND.I.NE.KHALF1)THEN
  194. YCOMP=CONJG(XX(NPOINT-I+2))
  195. PMODU2=PMODU2+XX(NPOINT-I+2)*YCOMP
  196. cbp: on pourrait aussi ecrire PMODU2=PMODU2+ABS(XX(NPOINT-I+2))**2
  197. ENDIF
  198. CPP +
  199. MLREE2.PROG(I-KDEBU+1)=COEF*PMODU2
  200. C IF(IIMPI.EQ.1) WRITE(IOIMP,*)'DSP(',I,')=',MLREE2.PROG(I)
  201. 20 CONTINUE
  202. C
  203. SEGDES MLREE1,MLREE2
  204. C
  205. C CREATION DE L'OBJET EVOLUTIO DSP
  206. C
  207. N=1
  208. SEGINI MEVOLL
  209. IPVO=MEVOLL
  210. TI(1:72)=TITREE
  211. IEVTEX=TI
  212. ITYEVO='REEL'
  213. SEGINI KEVOLL
  214. c KEVTEX=TI
  215. KEVTEX=KEVOL1.KEVTEX
  216. IEVOLL(1)=KEVOLL
  217. TYPX='LISTREEL'
  218. TYPY='LISTREEL'
  219. C
  220. IPROGX=MLREE1
  221. NOMEVX=MOTX(1:12)
  222. C
  223. IPROGY=MLREE2
  224. NOMEVY=MOTY(1:12)
  225. C
  226. NUMEVX=ICOUL
  227. NUMEVY='REEL'
  228. C
  229. SEGSUP MTRAV1,MTRAV2
  230. SEGDES KEVOL1
  231. SEGDES MEVOL1
  232. C
  233. SEGDES KEVOLL,MEVOLL
  234. CALL ECROBJ('EVOLUTIO',IPVO)
  235. 666 CONTINUE
  236. RETURN
  237. END
  238.  
  239.  
  240.  
  241.  
  242.  
  243.  
  244.  
  245.  
  246.  
  247.  
  248.  
  249.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales