anchor.diff:  
       use global_mod, only: dp, ip, force_double, charlen, wall_elevation, pi
       use domain_mod, only: domain_type, STG, UH, VH, ELV
  -    use read_raster_mod, only: read_gdal_raster
  +    use read_raster_mod, only: multi_raster_type
       use which_mod, only: which
       use ragged_array_mod, only: ragged_array_2d_ip_type
       use file_io_mod, only: read_character_file, read_csv_into_array
  
positive.code:         call read_csv_into_array(polygon_coords, polygon_filename)
          inside_point_counter = 0
          do j = 1, domain%nx(2)
  
              call points_in_poly(polygon_coords(1,:), polygon_coords(2,:), x(:,j), y(:,j), is_inside_poly)
              inside_point_counter = inside_point_counter + count(is_inside_poly)
  
              do i = 1, domain%nx(1) 
  
negative.code:         real(dp), parameter:: friction_road = 0.02_dp, friction_other = 0.04_dp
  
          ! things to help read polygons
          character(len=charlen), allocatable:: house_filenames(:), polygon_filename
          integer(ip):: house_file_unit, inside_point_counter, house_cell_count
          real(dp), allocatable:: polygon_coords(:,:)
          logical, allocatable:: is_inside_poly(:)
          
          ! For inflow boundary
          type(ragged_array_2d_ip_type), pointer :: rainfall_region_indices
  
  
          !
          ! Dry flow to start with. Later we clip stage to elevation.
  