Changeset 712 for trunk/LMDZ.MARS/libf


Ignore:
Timestamp:
Jun 29, 2012, 11:44:01 AM (12 years ago)
Author:
tnavarro
Message:

bug if icelocationmode=1 + sanity check for ice negative values

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/libf/phymars/surfini.F

    r697 r712  
    4545      INTEGER     imd,jmd
    4646      PARAMETER   (imd=360,jmd=180)
    47       REAL        zdata(imd*jmd)
     47      REAL        zdata(imd,jmd)
    4848      REAL        zelat,zelon
    49 ! on a data(i,j) = zdata((j-1)*jmd + i) avec i longitude de 1 a imd (E vers O) et j lat de 1 a jmd (N vers S)
    5049
    5150      INTEGER nb_ice(ngrid,2)              ! number of counts | detected ice for GCM grid
     
    190189      do ig=2,ngridmx-1
    191190     
    192         ! loop over the surface file grid
    193         do i=1,imd*jmd
    194          zelat = 90. - 90./real(jmd) - (i-1)/imd      ! surface file latitude
    195          zelon = mod(i-1,imd) - (180.-180./real(imd)) ! surface file longitude
    196                  
    197           if ((abs(lati(ig)*180./pi-zelat) .le. 90./real(jjm)) .and.
     191        ! loop over the surface file grid     
     192        do i=1,imd
     193          do j=1,jmd
     194            zelon = i - 180.
     195            zelat = 90. - j
     196            if ((abs(lati(ig)*180./pi-zelat) .le. 90./real(jjm)) .and.
    198197     &        (abs(long(ig)*180./pi-zelon) .le. 180./real(iim))) then
    199198              ! count all points in that GCM grid point
    200199              nb_ice(ig,1) = nb_ice(ig,1) + 1
    201               if (zdata(i) > min_icevalue)
     200              if (zdata(i,j) > min_icevalue)
    202201                 ! count all detected points in that GCM grid point
    203202     &           nb_ice(ig,2) = nb_ice(ig,2) + 1
    204           endif
    205          
    206         enddo ! of do i=1,imd*jmd
    207        
     203             endif
     204          enddo
     205        enddo  
     206
    208207        ! projection of nb_ice on GCM lat and lon axes
    209208        latice(1+(ig-2)/iim,:) =
     
    224223     
    225224     
    226   !    print*, 'latice TOT', latice(:,1)
    227   !    print*, 'latice FOUND', latice(:,2)
    228   !    print*, 'lonice TOT', lonice(:,1)
    229   !    print*, 'lonice FOUND', lonice(:,2)
    230      
    231   !    print*, 'lat ratio', int(real(latice(:,2))/real(latice(:,1))*iim)
    232   !    print*, 'lon ratio', int(real(lonice(:,2))/real(lonice(:,1))*jjm)
    233      
    234   !    print*,''
    235   !    print*,'sum lat', sum(latice(:,1)), sum(lonice(:,1))
    236   !    print*,'sum lon', sum(latice(:,2)), sum(lonice(:,2))
     225    print*, 'latice TOT', latice(:,1)
     226    print*, 'latice FOUND', latice(:,2)
     227    print*, 'lonice TOT', lonice(:,1)
     228    print*, 'lonice FOUND', lonice(:,2)
     229     
     230    print*, 'lat ratio', int(real(latice(:,2))/real(latice(:,1))*iim)
     231    print*, 'lon ratio', int(real(lonice(:,2))/real(lonice(:,1))*jjm)
     232     
     233    print*,''
     234    print*,'sum lat', sum(latice(:,1)), sum(lonice(:,1))
     235    print*,'sum lon', sum(latice(:,2)), sum(lonice(:,2))
    237236     
    238237   
     
    480479           if(noms(iq).eq."h2o_ice") then
    481480             do ig=1,ngrid
     481             
     482              if ((watercaptag(ig).eqv..false.)
     483     &     .and. (qsurf(ig,iq).lt.-frost_albedo_threshold)) then
     484              print*, ''
     485              print*, '!!! PROBLEM in SURFINI !!!!'
     486              print*, 'FOUND NEGATIVE SURFACE ICE VALUE WHERE
     487     & WATERCAPTAG IS FALSE'
     488              print*, ''
     489              print*, 'ig,qsurf,threshold' ,
     490     &         ig, qsurf(ig,iq), -frost_albedo_threshold
     491              print*, ''
     492              print*, '1) Check h2o_ice in startfi and ice
     493     & distribution in surfini'
     494              print*, '2) Use ini_h2osurf option in newstart'
     495              print*, ''
     496              CALL ABORT
     497            endif
     498             
    482499              if ((piceco2(ig) .eq. 0.).and.
    483500     &          (qsurf(ig,iq).gt.frost_albedo_threshold)) then
Note: See TracChangeset for help on using the changeset viewer.