source: trunk/northwind.pro @ 48

;
; NAME: northwind.pro
;
; PURPOSE: Check on in the intput field the level of coherence of winds nearby north pole,
;          and corrects it. Thi stask is performed before interpolation
;
; CATEGORY: Subroutine
;
; CALLING SEQUENCE: northwind
;
; INPUTS:
;
; KEYWORD PARAMETERS: None
;
; OUTPUTS:
;          zresul : the same wind field with its upward stripe
;          corrected in order to get coherence of winds on north polte
;
; COMMON BLOCKS:
;           common_interp.pro
;
PRO northwind, datglo,zdata_name,t

@common_interp

; Treatment of north pole stripe


; Get from global array, the last stripe and the last but one stripe

;    Initialising data arrays

     datglo_rect  = replicate(0.,2,jpiatm)
     datglo_rectb = replicate(0.,2,jpiatm)

; Put data of the last two stripes into arrays


     datglo_rect(0,*) = datglo(*,jpjatm-1,0)
     datglo_rect(1,*) = datglo(*,jpjatm-1,1)

     datglo_rectb(0,*)= datglo(*,jpjatm-2,0)
     datglo_rectb(1,*)= datglo(*,jpjatm-2,1)



; Converts data of these two lines into polar coordinates

     add_angle=0.

     if (zdata_name EQ data_u_name) then add_angle=!PI/2


     datglo_polar=CV_COORD(FROM_RECT=datglo_rect,/TO_POLAR)
     datglo_polar(0,*)=datglo_polar(0,*)+2*!PI*indgen(jpiatm)/jpiatm

     datglo_polarb=CV_COORD(FROM_RECT=datglo_rect,/TO_POLAR)
     datglo_polarb(0,*)=datglo_polar(0,*)+2*!PI*indgen(jpiatm)/jpiatm


; Coherence of the last line needs to be corrected

     print, 'CHECKING LAST LINE COHERENCE...'

     mean_mod=replicate(total(datglo_polar(1,*))/jpiatm,jpiatm)
     mean_arg=replicate((total(datglo_polar(0,*))/jpiatm+2*!PI*west_lon/360.+!PI*(jpiatm-1.)/float(jpiatm)) ,jpiatm)
     mean_argb=replicate((total(datglo_polarb(0,*))/jpiatm+2*!PI*west_lon/360.+!PI*(jpiatm-1.)/float(jpiatm)) ,jpiatm)



     ; Checks whether module is about the same on the North Pole line :

     test_mod=abs(mean_mod-datglo_polar(1,*))

     IF ( max(test_mod)/mean_mod(0) ) GT 0.01 THEN BEGIN

          print, 'MODULE IS NOT COHERENT AND WILL BE CORRECTED TO A MEAN VALUE OF :', mean_mod(0)

          datglo_polar(1,*)=mean_mod

     ENDIF


     ; Checks whether argument is coherent on the North Pole line :

     test_arg=abs(mean_arg-datglo_polar(0,*))

     IF ( max(test_arg)/mean_arg(0)) GT 0.01 THEN BEGIN

          print, 'ARGUMENT IS NOT COHERENT AND WILL BE CORRECTED TO A MEAN VALUE :', mean_argb(0)

          datglo_polar(0,*)=mean_argb-!PI/2-(2*!PI*indgen(jpiatm)/float(jpiatm)+replicate(2.*!PI*west_lon/360.,jpiatm))



     ; Cross Corellation Method for wind coherence on North Pole

     all_phases=indgen(jpiatm)

     Ux=mean_mod*cos(all_phases*2.*!PI/jpiatm)
     Uy=mean_mod*sin(all_phases*2.*!PI/jpiatm)

     phyx=C_CORRELATE(Ux,datglo(*,jpjatm-2,0),all_phases)

     phyy=C_CORRELATE(Uy,datglo(*,jpjatm-2,1),all_phases)


     deph_x=where(phyx EQ max(phyx))
     deph_y=where(phyy EQ max(phyy))

     deph_x=deph_x(0)
     deph_y=deph_y(0)




     Ux=shift(Ux,deph_x)
     Uy=shift(Uy,deph_y)

     datglo(*,jpjatm-1,0)=Ux
     datglo(*,jpjatm-1,1)=Uy


     ENDIF




end