Changeset 448 for trunk

Timestamp:

12/07/11 15:46:10 (13 years ago)

Author:

lelod

Message:

modifs extractamsua OK

Location:

trunk/src

Files:

• extract_amsua.pro (modified) (3 diffs)
• interpolswath.pro (modified) (3 diffs)

trunk/src/extract_amsua.pro

@@ -209 +209 @@
     openr, lun_a, list_filea, /get_lun
     filea=''
-    while (not eof(lun_a)) do begin
+    nlist=1
+    ilist=1
+   ; while (not eof(lun_a)) do begin
+    while ilist le nlist do begin
+      ilist=ilist+1
         onefile = ''
         readf, lun_a, filea
@@ -245 +249 @@
 ; latitude, et les Tbs des differents canaux (15 pour AMSUA)
 
-    amch1=REFORM(ama_scan.btemps[0,*,*]/100.)
-    amch2=REFORM(ama_scan.btemps[1,*,*]/100.)
-    amch3=REFORM(ama_scan.btemps[2,*,*]/100.)
-    amch4=REFORM(ama_scan.btemps[3,*,*]/100.)
-    amch5=REFORM(ama_scan.btemps[4,*,*]/100.)
-    amch6=REFORM(ama_scan.btemps[5,*,*]/100.)
-    amch7=REFORM(ama_scan.btemps[6,*,*]/100.)
-    amch8=REFORM(ama_scan.btemps[7,*,*]/100.)
-    amch9=REFORM(ama_scan.btemps[8,*,*]/100.)
-    amch10=REFORM(ama_scan.btemps[9,*,*]/100.)
-    amch11=REFORM(ama_scan.btemps[10,*,*]/100.)
-    amch12=REFORM(ama_scan.btemps[11,*,*]/100.)
-    amch13=REFORM(ama_scan.btemps[12,*,*]/100.)
-    amch14=REFORM(ama_scan.btemps[13,*,*]/100.)
-    amch15=REFORM(ama_scan.btemps[14,*,*]/100.)
+    ;amch1=REFORM(ama_scan.btemps[0,*,*]/100.)
+   ; amch2=REFORM(ama_scan.btemps[1,*,*]/100.)
+    ;amch3=REFORM(ama_scan.btemps[2,*,*]/100.)
+   ; amch4=REFORM(ama_scan.btemps[3,*,*]/100.)
+   ; amch5=REFORM(ama_scan.btemps[4,*,*]/100.)
+   ; amch6=REFORM(ama_scan.btemps[5,*,*]/100.)
+   ; amch7=REFORM(ama_scan.btemps[6,*,*]/100.)
+   ; amch8=REFORM(ama_scan.btemps[7,*,*]/100.)
+   ; amch9=REFORM(ama_scan.btemps[8,*,*]/100.)
+   ; amch10=REFORM(ama_scan.btemps[9,*,*]/100.)
+   ; amch11=REFORM(ama_scan.btemps[10,*,*]/100.)
+   ; amch12=REFORM(ama_scan.btemps[11,*,*]/100.)
+   ; amch13=REFORM(ama_scan.btemps[12,*,*]/100.)
+   ; amch14=REFORM(ama_scan.btemps[13,*,*]/100.)
+   ; amch15=REFORM(ama_scan.btemps[14,*,*]/100.)
+; extraction du canal traite
+    nocanal=5
+    amch=REFORM(ama_scan.btemps[nocanal,*,*]/100.)
     amzen=REFORM(ama_scan.angles[0,*,*]/100.)
     amalat=REFORM(ama_scan.latlon[0,*,*]/1.E4)
@@ -278 +285 @@
 ; appel a interpolswatw pour ajuster les pixels amsua sur une grille
 ; reguliere dans la fauchee
+; et selection de la zone conservee
 resol=1
-
+chint=fltarr(100)
+latint=fltarr(100)
+lonint=fltarr(100)
+fovint=fltarr(100)
+timeint=fltarr(100)
+time=fltarr(100)
 for i=0L,n_scan-1L do begin
-   tb=amch5[*,i]
+   tb=amch[*,i]
    lon=amalon[*,i]
    lat=amalat[*,i]
    interpolswath,tb,lat,lon,resol,nbgrid,tbgrid,latgrid,longrid
    if i eq 0L then begin
-     achint=fltarr(nbgrid,n_scan)
-     alatint=fltarr(nbgrid,n_scan)
-     alonint=fltarr(nbgrid,n_scan)
-     afovint=intarr(nbgrid,n_scan)
-     atimeint=fltarr(nbgrid,n_scan)
-     plot,lon,tb,psym=1
-     oplot,longrid,tbgrid,psym=4
+      fovgrid=indgen(nbgrid)+1
+      plot,lon,tb,psym=1
+      oplot,longrid,tbgrid,psym=4
   endif
-   achint[*,i]=tbgrid
-   alatint[*,i]=latgrid
-   alonint[*,i]=longrid
-   afovint[*,i]=indgen(nbgrid)+1
-   atimeint[*,i]=tt[i]
-endfor
-
-
-;-----------------
-; ici on selectionne la zone definie
-;-----------------
-    zone=where((alonint ge lon_min) and (alonint le lon_max) $
-             and (alatint ge lat_min) and (alatint le lat_max) )
-
-    nn=n_elements(zone)
-    if (nn GT 10) then begin
-        flagg=1
-        dam=fltarr(25,nn)
-        dam[0,*]=alonint[zone]
-        dam[1,*]=alatint[zone]
-       ; dam[2,*]=amazen[zone]
-        dam[3,*]=afovint[zone]
-        dam[4,*]=atimeint[zone]
-       ; dam[5,*]=amch1[zone]
-       ; dam[6,*]=amch2[zone]
-       ; dam[7,*]=amch3[zone]
-       ; dam[8,*]=amch4[zone]
-        dam[9,*]=achint[zone]
-       ; dam[10,*]=amch6[zone]
-      ;  dam[11,*]=amch7[zone]
-       ; dam[12,*]=amch8[zone]
-       ; dam[13,*]=amch9[zone]
-       ; dam[14,*]=amch10[zone]
-       ; dam[15,*]=amch11[zone]
-       ; dam[16,*]=amch12[zone]
-       ; dam[17,*]=amch13[zone]
-      ;  dam[18,*]=amch14[zone]
-      ;  dam[19,*]=amch15[zone]
- 
-
+  zone=where((longrid ge lon_min) and (longrid le lon_max) $
+             and (latgrid ge lat_min) and (latgrid le lat_max),npt)
+  if npt ne 0 then begin
+     time[zone]=tt[i]
+     chint=[chint,tbgrid[zone]]
+     latint=[latint,latgrid[zone]]
+     lonint=[lonint,longrid[zone]]
+     fovint=[fovint,fovgrid[zone]]
+     timeint=[timeint,time[zone]]
+  endif
+
+ endfor
+
+ ; partie a remplacer par module ecriture standardise FP
+;
       ; replace missing values -9999.99 by NaN
-        index_missing = where(dam EQ -9999.99, nb_missing)
-        if (nb_missing NE 0) THEN BEGIN
-            dam[index_missing]=!values.f_nan
-        ENDIF
-
-      ; ecriture des donnees ligne par ligne dans le fichier resultant
-      ;
-      ; mm : mois, a recuperer dans les parametres
-      ; nosat: 15, 16 ou 17
-     ;   print,'orbite',desc
-      for i=0L, nn-1L do begin
-          
+  index_missing = where(chint EQ -9999.99, nb_missing)
+  if (nb_missing NE 0) THEN BEGIN
+     chint[index_missing]=!values.f_nan
+  endif
+; param landseamask mis fictivement pour unification du format
+        landseamask=!values.f_nan
+nn=max(size(chint))
+print,'nb points dans la zone',nn
+
+      for i=0L, nn-1L do begin     
 ;format: nosat,mois, day, lon, lat, zen,fov, time, ama_ch1:15,amb_ch1:5
- ; print, 'mois::::', mm,'-----', nab-i
-              printf,lun_obs, mm, jour_ama,desc, dam[4,i], nosat, dam[3,i], dam[0,i], dam[1,i],$
-                dam[9,i],format='(3(i3,1x),f8.4,2i3,3(1x,f8.3))'
+; print, 'mois::::', mm,'-----', nab-i
+         printf,lun_obs, mm, jour_ama, desc, nosat, fovint[i], landseamask, lonint[i], latint[i], timeint[i],$
+                chint[i],format='(6(i3,1x),4(f8.4,1x))'
      endfor
       
-  endif else begin
-      flagg=0
-      print, 'pas de points dans la zone'
-  endelse
 
 

trunk/src/interpolswath.pro

@@ -56 +56 @@
 
 for i=nnadir+1,na-1 do begin
-dist=rterre*acos(sin(lat(i-1))*sin(lat(i)) + cos(lat(i-1))*cos(lat(i))*cos(lon(i)-lon(i-1)))
-fova[i]=fova[i-1]+dist
+   dist=rterre*acos(sin(lat(i-1))*sin(lat(i)) + cos(lat(i-1))*cos(lat(i))*cos(lon(i)-lon(i-1)))
+   fova[i]=fova[i-1]+dist
 endfor
 for i=nnadir-2,0,-1 do begin
-dist=rterre*acos(sin(lat(i))*sin(lat(i+1)) + cos(lat(i))*cos(lat(i+1))*cos(lon(i)-lon(i+1)))
+   dist=rterre*acos(sin(lat(i))*sin(lat(i+1)) + cos(lat(i))*cos(lat(i+1))*cos(lon(i)-lon(i+1)))
 ;print,i,dist
-fova[i]=fova[i+1]-dist
+   fova[i]=fova[i+1]-dist
 endfor
 
 print,fova
-;if -fova[0] ne fova[na-1] then print, 'pb symetrie fauchee'
 
 distmax=-min(fova)+max(fova)+max(pixa)
@@ -79 +78 @@
 grid2=indgen(ndemigrid)*delta+fova[nnadir]
 grid1=fltarr(ndemigrid)
-for i=0,ndemigrid-1 do grid1[i]=-grid2[ndemigrid-1-i]
+for i=0,ndemigrid-1 do begin
+   grid1[i]=-grid2[ndemigrid-1-i]
+endfor
 ;alternative equivalente)
 ;for i=ndemigrid-1,0,-1 do grid1[i]=-(ndemigrid-i-1)*delta-fova[nnadir]
@@ -87 +88 @@
 tbint=fltarr(nbgrid)
 for i=0,na-1 do begin
-ind=where(abs(grid-fova[i]) le pixatot[i]/2,nii)
-;print,i,nii
-tbint[ind]=tb[i]
+   ind=where(abs(grid-fova[i]) le pixatot[i]/2,nii)
+   print,i,nii
+   if nii ne 0 then begin
+      tbint[ind]=tb[i]
+   endif
 endfor
 latgrid=interpol(latlu,fova,grid)