anchor.diff:      real,parameter :: freezing_threshold=273.15          ! [K]
       real,parameter :: snow_fall_rate=1.5                 ! [m/s]   for a water vapor scale height of 3750m corresponds to tau_f = 2500
       real,parameter :: rain_fall_rate=10.0                ! [m/s]   for a water vapor scale height of 3750m corresponds to tau_f = 375
  -    real,parameter :: snow_cloud_init=0.00001            ! [kg/kg] cloud ice content before snow will start to form
  -    real,parameter :: rain_cloud_init=0.00001            ! [kg/kg] cloud water content before rain will start to form
  +    real,parameter :: snow_cloud_init=0.0001            ! [kg/kg] cloud ice content before snow will start to form
  +    real,parameter :: rain_cloud_init=0.0001            ! [kg/kg] cloud water content before rain will start to form
       
       
   !   these are recalculated every call because they are a function of dt
  
positive.code:             if (qv<qvsat) then
                  if (qr>SMALL_VALUE) then
                      ! evaporate rain
                      call phase_change(p,t,qr,qvsat,qv,L_evap,1.0) !cloud2rain)
                  endif
                  if (qs>SMALL_VALUE) then
                      ! sublimate snow
                      call phase_change(p,t,qs,qvsat,qv,L_subl,1.0) !cloud2snow)
                  endif
              endif
          endif
  
negative.code:                 th(i,:,j)=t/pii(i,:,j)
                  
              enddo
  !             where(qs(:,:,j)<0) qs(:,:,j)=0
  !             where(qr(:,:,j)<0) qr(:,:,j)=0
  !             where(qc(:,:,j)<0) qc(:,:,j)=0
  !             where(qv(:,:,j)<0) qv(:,:,j)=0
          enddo
          !$omp end do
          deallocate(t)
  