- Timestamp:
- Sep 18, 2018 9:53:18 AM (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/UTIL/agent_preprocessing/agent_preprocessing.f90
r3216 r3259 25 25 ! -----------------! 26 26 ! $Id$ 27 ! Removed unused variables and fixed real to real comparisons 28 ! 29 ! 3216 2018-08-29 10:22:12Z sward 27 30 ! Bubfix for gfortran: reordering of type definitions 28 31 ! … … 248 251 249 252 INTEGER(iwp) :: i !< running index along x-direction 250 INTEGER(iwp) :: ii !< running index for IO blocks251 253 INTEGER(iwp) :: k_head !< minimum k index for agents to walk underneath overhanging buildings 252 254 INTEGER(iwp) :: id_topo !< NetCDF id of topograhy input file 253 255 INTEGER(iwp) :: j !< running index along y-direction 254 INTEGER(iwp) :: k !< running index along z-direction255 256 INTEGER(iwp) :: num_vars !< number of variables in netcdf input file 256 INTEGER(iwp) :: skip_n_rows !< counting variable to skip rows while reading topography file257 257 258 258 LOGICAL :: netcdf_flag = .FALSE. !< indicates whether netcdf file is used for input 259 259 LOGICAL :: lod_flag = .FALSE. !< true if 3d building data is used 260 260 LOGICAL :: topo_file_flag = .FALSE. !< true if 3d building data is used 261 262 REAL(wp) :: dum !< dummy variable to skip columns while reading topography file263 261 264 262 WRITE(*,'((X,A,/))') 'Looking for topography/building information' … … 514 512 CHARACTER (LEN=*), INTENT(IN) :: filename !< filename 515 513 INTEGER(iwp), INTENT(INOUT) :: id !< file id 516 LOGICAL :: file_open = .FALSE.517 514 518 515 nc_stat = NF90_OPEN( filename, NF90_NOWRITE, id ) … … 803 800 INTEGER(iwp), INTENT(IN) :: n3 !< number of data-points along 3rd dimension 804 801 INTEGER(iwp), INTENT(IN) :: n4 !< number of data-points along 4th dimension 805 806 INTEGER(iwp), DIMENSION(3) :: id_dim807 802 808 803 REAL(wp), DIMENSION(:,:), INTENT(INOUT) :: var !< variable to be read … … 1242 1237 IMPLICIT NONE 1243 1238 1244 CHARACTER(LEN=20) :: FMT1245 1239 CHARACTER(LEN=200) :: dirname 1246 1240 CHARACTER(LEN=200) :: rundir … … 1355 1349 !-- If tolerance_dp was not set, put in default values 1356 1350 DO i = 0, 2 1357 IF ( tolerance_dp(i) == 999999.0_wp ) THEN1351 IF ( tolerance_dp(i) > 999998.0_wp ) THEN 1358 1352 tolerance_dp(i) = SQRT(dx*dy)*1.41/(2**i) 1359 1353 ELSE … … 1368 1362 ! 1369 1363 !-- Set null_vertex 1370 CALL set_vertex(null_vertex,0 _iwp,0.0_wp,0.0_wp)1364 CALL set_vertex(null_vertex,0.0_wp,0.0_wp) 1371 1365 ! 1372 1366 !-- Some initializations … … 1571 1565 IMPLICIT NONE 1572 1566 1573 INTEGER(iwp) :: il !< local counter1574 1567 INTEGER(iwp) :: p_id !< current polygon_id 1575 1568 ! … … 1588 1581 polygons(p_id)%nov = polygons(p_id)%nov + 1 1589 1582 nov = polygons(p_id)%nov 1590 CALL set_vertex(dummy_vertex, p_id,i*dx, j*dy)1583 CALL set_vertex(dummy_vertex, i*dx, j*dy) 1591 1584 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1592 1585 ENDIF … … 1599 1592 polygons(p_id)%nov = polygons(p_id)%nov + 1 1600 1593 nov = polygons(p_id)%nov 1601 CALL set_vertex(dummy_vertex, p_id,i*dx, (j+1)*dy)1594 CALL set_vertex(dummy_vertex, i*dx, (j+1)*dy) 1602 1595 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1603 1596 ENDIF … … 1610 1603 polygons(p_id)%nov = polygons(p_id)%nov + 1 1611 1604 nov = polygons(p_id)%nov 1612 CALL set_vertex(dummy_vertex, p_id,(i+1)*dx, (j+1)*dy)1605 CALL set_vertex(dummy_vertex, (i+1)*dx, (j+1)*dy) 1613 1606 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1614 1607 ENDIF … … 1621 1614 polygons(p_id)%nov = polygons(p_id)%nov + 1 1622 1615 nov = polygons(p_id)%nov 1623 CALL set_vertex(dummy_vertex, p_id,(i+1)*dx, j*dy)1616 CALL set_vertex(dummy_vertex, (i+1)*dx, j*dy) 1624 1617 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1625 1618 ENDIF … … 1636 1629 polygons(p_id)%nov = polygons(p_id)%nov + 1 1637 1630 nov = polygons(p_id)%nov 1638 CALL set_vertex(dummy_vertex, p_id,i*dx, j*dy)1631 CALL set_vertex(dummy_vertex, i*dx, j*dy) 1639 1632 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1640 1633 ENDIF … … 1648 1641 polygons(p_id)%nov = polygons(p_id)%nov + 1 1649 1642 nov = polygons(p_id)%nov 1650 CALL set_vertex(dummy_vertex, p_id,i*dx, (j+1)*dy)1643 CALL set_vertex(dummy_vertex, i*dx, (j+1)*dy) 1651 1644 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1652 1645 ENDIF … … 1660 1653 polygons(p_id)%nov = polygons(p_id)%nov + 1 1661 1654 nov = polygons(p_id)%nov 1662 CALL set_vertex(dummy_vertex, p_id,(i+1)*dx, (j+1)*dy)1655 CALL set_vertex(dummy_vertex, (i+1)*dx, (j+1)*dy) 1663 1656 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1664 1657 ENDIF … … 1672 1665 polygons(p_id)%nov = polygons(p_id)%nov + 1 1673 1666 nov = polygons(p_id)%nov 1674 CALL set_vertex(dummy_vertex, p_id,(i+1)*dx, j*dy)1667 CALL set_vertex(dummy_vertex, (i+1)*dx, j*dy) 1675 1668 CALL add_vertex_to_polygon(dummy_vertex, p_id, nov) 1676 1669 ENDIF … … 1686 1679 !> Initializes a vertex 1687 1680 !------------------------------------------------------------------------------! 1688 SUBROUTINE set_vertex (in_vertex, p_id, x, y) 1689 1690 IMPLICIT NONE 1691 1692 INTEGER(iwp) :: p_id !< polygon ID 1681 SUBROUTINE set_vertex (in_vertex, x, y) 1682 1683 IMPLICIT NONE 1693 1684 1694 1685 REAL(wp) :: x !< x-coordinate of vertex position … … 1798 1789 ENDIF 1799 1790 DO il = 1, nov 1800 IF ( sorted_p(il)%x == m_x ) THEN1791 IF ( ABS( sorted_p(il)%x - m_x ) < .01 * dx ) THEN 1801 1792 counter = counter + 1 1802 1793 dummy_vertex = sorted_p(il) … … 1811 1802 ENDIF 1812 1803 DO il = 1, counter 1813 IF ( sorted_p(il)%y == m_y ) THEN1804 IF ( ABS(sorted_p(il)%y - m_y) < .01 * dy ) THEN 1814 1805 dummy_vertex = sorted_p(il) 1815 1806 sorted_p(il) = sorted_p(1) … … 1860 1851 IF ( nosv < nov ) THEN 1861 1852 DO il = nosv+1, nov 1862 IF ( current_v%x == sorted_p(il)%x .OR.&1863 current_v%y == sorted_p(il)%y)&1853 IF ( ABS( current_v%x - sorted_p(il)%x ) < .01 * dx .OR. & 1854 ABS( current_v%y - sorted_p(il)%y ) < .01 * dy ) & 1864 1855 THEN 1865 1856 ! … … 1888 1879 !-- all tests fail, is the candidate one of a maximum of two possible 1889 1880 !-- neighbors. 1881 id_neighbor = -999 1890 1882 id_neighbor1 = -999 1891 1883 id_neighbor2 = -999 … … 1907 1899 (candidates(il)%y - current_v%y)**2 ) 1908 1900 IF ( nosv > 1 ) THEN 1909 IF ( dist < .9*dx .OR. & 1910 (sorted_p(nosv-1)%x == candidates(il)%x .AND. & 1911 sorted_p(nosv-1)%y == candidates(il)%y) ) & 1901 IF ( dist < .9 * dx .OR. & 1902 ( ( ABS( sorted_p(nosv-1)%x & 1903 - candidates(il)%x ) < .01 * dx ) .AND. & 1904 ( ABS( sorted_p(nosv-1)%y & 1905 - candidates(il)%y ) < .01 * dy ) ) ) & 1912 1906 THEN 1913 1907 CYCLE … … 1919 1913 !-- (4 cases) 1920 1914 !-- First: for vertical connection 1921 IF ( candidates(il)%x == current_v%x ) THEN1915 IF ( ABS( candidates(il)%x - current_v%x ) < .01 * dx ) THEN 1922 1916 xs = NINT(current_v%x*ddx)-1 1923 1917 xe = xs + 1 … … 1947 1941 ! 1948 1942 !-- Horizontal connection 1949 ELSEIF ( candidates(il)%y == current_v%y ) THEN1943 ELSEIF ( ABS( candidates(il)%y - current_v%y ) < .01 * dy ) THEN 1950 1944 1951 1945 ys = NINT(current_v%y*ddy)-1 … … 1982 1976 id_neighbor1 = candidate_id(il) 1983 1977 ELSEIF ( id_neighbor1 /= -999 .AND. & 1984 ( sorted_p(id_neighbor1)%x /= candidates(il)%x ) .OR. & 1985 ( sorted_p(id_neighbor1)%y /= candidates(il)%y ) ) & 1978 ( ( ABS( sorted_p(id_neighbor1)%x - candidates(il)%x ) & 1979 > .01 * dx ) .OR. & 1980 ( ABS( sorted_p(id_neighbor1)%y - candidates(il)%y ) & 1981 > .01 * dy ) ) ) & 1986 1982 THEN 1987 1983 id_neighbor2 = candidate_id(il) … … 2321 2317 INTEGER(iwp) :: vid_t !< vertex id of tested mesh point 2322 2318 INTEGER(iwp) :: vl !< vertex counter 2323 2324 REAL(wp) :: left !< counter: current mesh point 2319 2325 2320 REAL(wp) :: v1x !< x-coordinate of test vertex 1 for intersection test 2326 2321 REAL(wp) :: v1y !< y-coordinate of test vertex 1 for intersection test … … 2678 2673 IMPLICIT NONE 2679 2674 2680 INTEGER(iwp) :: il !< Local counter2681 2675 INTEGER(iwp) :: noc !< Number of connections 2682 2676
Note: See TracChangeset
for help on using the changeset viewer.