anchor.diff:                print*, "    Using a sleve_n of ", options%parameters%sleve_n
                 ! print*, ""
                 write(*,*) "    Smooth height (model top) is ", smooth_height, "m.a.s.l"
  +              write(*,*) "    invertibility parameter gamma is: ", gamma
  +              if(gamma <= 0) print*, " CAUTION: coordinate transformation is not invertible (gamma <= 0 ) !!! reduce decay rate(s)!"
                 ! write(*,*) "  mean terrain ", sum(terrain) / size(terrain)
                 ! write(*,*) "  sum(dz) ", sum(dz(1:max_level))
                 ! write(*,*) "  sum(dz_scl) ", sum(dz_scl(1:max_level))
  
positive.code:             ! H        =  sum(dz_scl(1:max_level))  ! should also lead to smooth_height, but more error proof?
  
  
              ! Decay Rate for Large-Scale Topography: svc1 = 10000.0000  COSMO1 operational setting (but model top is at ~22000 masl)
              ! Decay Rate for Small-Scale Topography: svc2 =  3300.0000
              if ((this_image()==1)) then
  
negative.code:                   jms => this%jms,      jme => this%jme,                        &
                    kms => this%kms,      kme => this%kme,                        &
                    terrain               => this%terrain%data_2d,                &
                    global_terrain        => this%global_terrain,                 &
                    terrain_u             => this%terrain_u,                      &
                    terrain_v             => this%terrain_v,                      &
                    forcing_terrain       => this%forcing_terrain%data_2d,        &
  