- Timestamp:
- 05/02/06 14:11:53 (18 years ago)
- Location:
- trunk
- Files:
-
- 1 deleted
- 3 copied
-
STATISTICS (deleted)
-
ToBeReviewed/STATISTICS/a_correlate2d.pro (copied) (copied from trunk/STATISTICS/a_correlate2d.pro)
-
ToBeReviewed/STATISTICS/a_timecorrelate.pro (copied) (copied from trunk/STATISTICS/a_timecorrelate.pro)
-
ToBeReviewed/STATISTICS/c_timecorrelate.pro (copied) (copied from trunk/STATISTICS/c_timecorrelate.pro) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/ToBeReviewed/STATISTICS/c_timecorrelate.pro
r19 r21 97 97 ; 98 98 ; MODIFICATION HISTORY: 99 ; 01/03/2000 Sebastien Masson (smasson@lodyc.jussieu.fr)99 ; - 01/03/2000 Sebastien Masson (smasson@lodyc.jussieu.fr) 100 100 ; Based on the C_CORRELATE procedure of IDL 101 ; - August 2003 Sebastien Masson 102 ; update according to the update made in C_CORRELATE by 103 ; W. Biagiotti and available in IDL 5.5 101 104 ;- 102 105 103 FUNCTION TimeCross_Cov, X, Y, M, nT, Double = Double, ZERO2NAN = zero2nan 104 ; 105 if double then one = 1.0d ELSE one = 1.0 106 ; Sample cross covariance function. 107 TimeDim = size(X, /n_dimensions) 108 Xmean = TOTAL(X, TimeDim, Double = Double) / nT 109 Xmean = Xmean[*]#replicate(one, nT - M) 110 Ymean = TOTAL(Y, TimeDim, Double = Double) / nT 111 Ymean = Ymean[*]#replicate(one, nT - M) 112 ; 113 case TimeDim of 114 1:res = TOTAL((X[0:nT - M - 1L] - Xmean) * (Y[M:nT - 1L] - Ymean) $ 106 FUNCTION TimeCross_Cov, Xd, Yd, M, nT, Ndim, Double = Double, ZERO2NAN = zero2nan 107 ;Sample cross covariance function. 108 109 compile_opt hidden 110 ; 111 case Ndim OF 112 1:res = TOTAL(Xd[0:nT - M - 1L] * Yd[M:nT - 1L] $ 115 113 , Double = Double) 116 2:res = TOTAL( (X[*, 0:nT - M - 1L] - Xmean) * (Y[*, M:nT - 1L] - Ymean)$117 , TimeDim, Double = Double)118 3:res = TOTAL( (X[*, *, 0:nT - M - 1L] - Xmean) * (Y[*, *, M:nT - 1L] - Ymean)$119 , TimeDim, Double = Double)120 4:res = TOTAL( (X[*, *, *, 0:nT - M - 1L] - Xmean) * (Y[*, *, *, M:nT - 1L] - Ymean)$121 , TimeDim, Double = Double)114 2:res = TOTAL(Xd[*, 0:nT - M - 1L] * Yd[*, M:nT - 1L] $ 115 , Ndim, Double = Double) 116 3:res = TOTAL(Xd[*, *, 0:nT - M - 1L] * Yd[*, *, M:nT - 1L] $ 117 , Ndim, Double = Double) 118 4:res = TOTAL(Xd[*, *, *, 0:nT - M - 1L] * Yd[*, *, *, M:nT - 1L] $ 119 , Ndim, Double = Double) 122 120 ENDCASE 123 121 if keyword_set(zero2nan) then begin … … 157 155 nLag = N_ELEMENTS(Lag) 158 156 157 ;Deviations 158 if double then one = 1.0d ELSE one = 1.0 159 Ndim = size(X, /n_dimensions) 160 Xd = TOTAL(X, Ndim, Double = Double) / nT 161 Xd = X - Xd[*]#replicate(one, nT) 162 Yd = TOTAL(Y, Ndim, Double = Double) / nT 163 Yd = Y - Yd[*]#replicate(one, nT) 164 159 165 if nLag eq 1 then Lag = [Lag] ;Create a 1-element vector. 160 166 161 167 case NDim of 162 1:if Double eq 0 then Correl = FLTARR(nLag) else Correl= DBLARR(nLag)163 2:if Double eq 0 then Correl = FLTARR(Xsize[1], nLag) else Correl= DBLARR(Xsize[1], nLag)164 3:if Double eq 0 then Correl= FLTARR(Xsize[1], Xsize[2], nLag) $165 else Correl= DBLARR(Xsize[1], Xsize[2], nLag)166 4:if Double eq 0 then Correl= FLTARR(Xsize[1], Xsize[2], Xsize[3], nLag) $167 else Correl= DBLARR(Xsize[1], Xsize[2], Xsize[3], nLag)168 1:if Double eq 0 then Cross = FLTARR(nLag) else Cross = DBLARR(nLag) 169 2:if Double eq 0 then Cross = FLTARR(Xsize[1], nLag) else Cross = DBLARR(Xsize[1], nLag) 170 3:if Double eq 0 then Cross = FLTARR(Xsize[1], Xsize[2], nLag) $ 171 else Cross = DBLARR(Xsize[1], Xsize[2], nLag) 172 4:if Double eq 0 then Cross = FLTARR(Xsize[1], Xsize[2], Xsize[3], nLag) $ 173 else Cross = DBLARR(Xsize[1], Xsize[2], Xsize[3], nLag) 168 174 endcase 169 175 170 if KEYWORD_SET(Covariance) eq 0 then begin ;Compute Cross Correlation.176 if KEYWORD_SET(Covariance) eq 0 then begin ;Compute Cross Crossation. 171 177 for k = 0, nLag-1 do begin 172 178 if Lag[k] ge 0 then BEGIN 173 179 case NDim of 174 1:Correl[k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / $ 175 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 176 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 177 2:Correl[*, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / $ 178 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 179 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 180 3:Correl[*, *, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / $ 181 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 182 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 183 4:Correl[*, *, *, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / $ 184 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 185 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 186 endcase 180 1: Cross[k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) 181 2: Cross[*, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) 182 3: Cross[*, *, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) 183 4: Cross[*, *, *, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) 184 endcase 187 185 ENDIF else BEGIN 188 186 case NDim of 189 1:Correl[k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / $ 190 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 191 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 192 2:Correl[*, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / $ 193 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 194 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 195 3:Correl[*, *, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / $ 196 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 197 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 198 4:Correl[*, *, *, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / $ 199 sqrt(TimeCross_Cov(X,X, 0L, nT, Double = Double, /zero2nan) * $ 200 TimeCross_Cov(Y,Y, 0L, nT, Double = Double, /zero2nan)) 201 endcase 187 1: Cross[k] = TimeCross_Cov(Yd, Xd, ABS(Lag[k]), nT, Ndim, Double = Double) 188 2: Cross[*, k] = TimeCross_Cov(Yd, Xd, ABS(Lag[k]), nT, Ndim, Double = Double) 189 3: Cross[*, *, k] = TimeCross_Cov(Yd, Xd, ABS(Lag[k]), nT, Ndim, Double = Double) 190 4: Cross[*, *, *, k] = TimeCross_Cov(Yd, Xd, ABS(Lag[k]), nT, Ndim, Double = Double) 191 endcase 202 192 ENDELSE 203 endfor 193 ENDFOR 194 div = sqrt(TimeCross_Cov(Xd, Xd, 0L, nT, Ndim, Double = Double, /zero2nan) * $ 195 TimeCross_Cov(Yd, Yd, 0L, nT, Ndim, Double = Double, /zero2nan)) 196 Cross = temporary(Cross)/((temporary(div))[*]#replicate(one, nLag)) 204 197 endif else begin ;Compute Cross Covariance. 205 198 for k = 0, nLag-1 do begin 206 199 if Lag[k] ge 0 then BEGIN 207 200 case NDim of 208 1: Correl[k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / nT209 2: Correl[*, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / nT210 3: Correl[*, *, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / nT211 4: Correl[*, *, *, k] = TimeCross_Cov(X, Y, Lag[k], nT, Double = Double) / nT201 1: Cross[k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) / nT 202 2: Cross[*, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) / nT 203 3: Cross[*, *, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) / nT 204 4: Cross[*, *, *, k] = TimeCross_Cov(Xd, Yd, Lag[k], nT, Ndim, Double = Double) / nT 212 205 ENDCASE 213 206 ENDIF else BEGIN 214 207 case NDim of 215 1: Correl[k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / nT216 2: Correl[*, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / nT217 3: Correl[*, *, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / nT218 4: Correl[*, *, *, k] = TimeCross_Cov(y, x, ABS(Lag[k]), nT, Double = Double) / nT208 1: Cross[k] = TimeCross_Cov(yd, xd, ABS(Lag[k]), nT, Ndim, Double = Double) / nT 209 2: Cross[*, k] = TimeCross_Cov(yd, xd, ABS(Lag[k]), nT, Ndim, Double = Double) / nT 210 3: Cross[*, *, k] = TimeCross_Cov(yd, xd, ABS(Lag[k]), nT, Ndim, Double = Double) / nT 211 4: Cross[*, *, *, k] = TimeCross_Cov(yd, xd, ABS(Lag[k]), nT, Ndim, Double = Double) / nT 219 212 ENDCASE 220 213 ENDELSE … … 222 215 endelse 223 216 224 if Double eq 0 then RETURN, FLOAT(Correl) else $ 225 RETURN, Correl 217 if Double eq 0 then RETURN, FLOAT(Cross) else RETURN, Cross 226 218 227 219 END
Note: See TracChangeset
for help on using the changeset viewer.
