Skip to content
Snippets Groups Projects
Commit acd9d06b authored by Anthony's avatar Anthony
Browse files

Improvement of Linkup

parent 24a74187
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
F90subroutines/routing_reg.f90
+ 67
− 58
View file @ acd9d06b
...@@ -2054,28 +2054,30 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area, ...@@ -2054,28 +2054,30 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area,
IF ( outflow_basin(sp,sb) == bop ) THEN IF ( outflow_basin(sp,sb) == bop ) THEN
solved(sp,1) = solved(sp,1) + 1 solved(sp,1) = solved(sp,1) + 1
ENDIF ENDIF
WRITE(numout,*) sp,sb,"flow in the same grid, basin:",bop
! !
ELSE IF ( outflow_grid(sp,sb) .EQ. -3 ) THEN ELSE IF ( outflow_grid(sp,sb) .EQ. -3 ) THEN
! Return flow ! Return flow
! Nothing to do but just remember it is done. ! Nothing to do but just remember it is done.
solved(sp,1) = solved(sp,1) + 1 solved(sp,1) = solved(sp,1) + 1
WRITE(numout,*) sp,sb,"is a return flow"
ELSE IF ( outflow_grid(sp,sb) .EQ. -2 ) THEN ELSE IF ( outflow_grid(sp,sb) .EQ. -2 ) THEN
! Coastal flow ! Coastal flow
! Nothing to do but just remember it is done. ! Nothing to do but just remember it is done.
solved(sp,1) = solved(sp,1) + 1 solved(sp,1) = solved(sp,1) + 1
WRITE(numout,*) sp,sb,"is a coastal flow"
ELSE IF ( outflow_grid(sp,sb) .EQ. -1 ) THEN ELSE IF ( outflow_grid(sp,sb) .EQ. -1 ) THEN
! River flow ! River flow
! Nothing to do but just remember it is done. ! Nothing to do but just remember it is done.
solved(sp,1) = solved(sp,1) + 1 solved(sp,1) = solved(sp,1) + 1
WRITE(numout,*) sp,sb,"is a river outflow"
ENDIF ENDIF
WRITE(numout,*) sp,sb,"is an outflow ", outflow_grid(sp,sb)
ELSE IF ( outflow_grid(sp,sb) .GT. 0 ) THEN ELSE IF ( outflow_grid(sp,sb) .GT. 0 ) THEN
found = 0 found = 0
! !
! Deal with the first one as usual ! Deal with the first one as usual
! !
inp = outflow_grid(sp,sb) inp = outflow_grid(sp,sb)
IF ( inp == 0.) WRITE(numout,*) "ERROR 1"
! !
CALL routing_reg_bestsubbasin(sp, sb, basin_id(sp,sb), basin_hierarchy(sp,sb), basin_fac(sp,sb), & CALL routing_reg_bestsubbasin(sp, sb, basin_id(sp,sb), basin_hierarchy(sp,sb), basin_fac(sp,sb), &
& basin_flowdir(sp,sb), invented_basins, & & basin_flowdir(sp,sb), invented_basins, &
...@@ -2223,27 +2225,10 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area, ...@@ -2223,27 +2225,10 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area,
! ! loc(NbNeighb) ! ! loc(NbNeighb)
ENDIF ENDIF
! !
! If not found, and at least one of the neighbour is ocean -> coastal outflow ! Looking for a solution in the grid -> HTU with a similar hierarchy or lowest hierarchy
IF ( found == 0) THEN !
WRITE (numout,*) "Is there a coastal outflow ? "
nbocean = 0
DO nb=1,NbNeighb
IF ( neighbours_g(sp,nb) .LE. -1 ) nbocean = nbocean +1
ENDDO
IF ( nbocean .LE. 1 ) THEN
outflow_grid(sp,sb) = -2
outflow_basin(sp,sb) = undef_int
found = 1
WRITE (numout,*) "There is more than one neighbour that is coastal outflow"
WRITE (numout,*) "Changed in coastal outflow"
solved(sp,3) = solved(sp,3)+1
ENDIF
ENDIF
! Equivalent to linkup 5
IF ( found == 0 ) THEN IF ( found == 0 ) THEN
WRITE (numout,*) "Eq. Linkup 5" WRITE (numout,*) "Looking for a solution in the grid"
DO sba=1,basin_count(sp) DO sba=1,basin_count(sp)
IF ( sba .NE. sb ) THEN IF ( sba .NE. sb ) THEN
IF ( basin_id(sp,sb) == basin_id(sp,sba) ) THEN IF ( basin_id(sp,sb) == basin_id(sp,sba) ) THEN
...@@ -2285,47 +2270,71 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area, ...@@ -2285,47 +2270,71 @@ SUBROUTINE routing_reg_linkup(nbpt, neighbours, nwbas, basin_count, basin_area,
sbint = minhloc(basin_id(sp,sb),2) sbint = minhloc(basin_id(sp,sb),2)
ENDIF ENDIF
ENDIF ENDIF
!
! Trung: Be careful when changing value of basin_hierarchy here ! IF (( sbint < undef_int ) .AND. (sbint .GT. 0) .AND. (sbint .NE. sb)) THEN
!
IF (( sbint < undef_int ) .AND. (sbint .GT. 0) ) THEN
basin_hierarchy(sp,sb) = basin_hierarchy(sp,sbint) basin_hierarchy(sp,sb) = basin_hierarchy(sp,sbint)
CALL routing_updateflow(sp, sb, sp, sbint, nbpt, nwbas, nbxmax*8, outflow_grid, outflow_basin, & CALL routing_updateflow(sp, sb, sp, sbint, nbpt, nwbas, nbxmax*8, outflow_grid, outflow_basin, &
& inflow_number, inflow_grid, inflow_basin) & inflow_number, inflow_grid,inflow_basin)
found = 1
WRITE (numout,*) "Lowest basin hierarchy in the grid file" WRITE (numout,*) "Lowest basin hierarchy in the grid file"
ELSE
!
! Trung: Well, if you come to this step, it should be a special
! (or curious, I said) case. For example, when processing
! Sulina branch of river Danube with the forcing data
! from E2OFD (/bdd/ORCHIDEE_Forcing/BC/OOL/OL2/E2OFD/)
! in resolution of 0.25 degree. The grid box [45.125N, 29.375E]
! is a ocean point in the middle of 7 land points ("a
! hole"). But HydroSHEDS shows that Sulina branch flows
! through this grid cell (from [45.125N, 29.125E] to [45.125N, 29.625E]).
! Trung: If the "JUMP" does not work, it means a worse case.
! Or it means a worse case.
! The outlet of this river falls in the ocean point of
! forcing data. For example: with 0.25 degree forcing
! data of E2OFD, the grid point [ 45.625N ; 13.875E ]
! lost the outlet of basin ID 69666.
! We need to use this point (which should be the last point
! of this river with shortest hierarchy and highest flow accumulation)
! to become new outlet.
!
! This is not good solution because there can be few
! outlet for one basin ID.
! But I need the model works now !!! => so, come back
! here later !
!
WRITE (numout,*) "Linkup : Made a NEW OUTLET at sp & sb: ",sp,sb
! Coastal flow or river flow is both ok here
outflow_grid(sp,sb) = -2
basin_hierarchy(sp,sb) = 0.00
ENDIF ENDIF
IF ( outflow_basin(sp,sb) < undef_int ) solved(sp,4) = solved(sp,4)+1
ENDIF ENDIF
!
! Looking if the minimun hierarchy is in a neighbour grid point
IF (found .EQ. 0) THEN
dop = undef_int
bop = undef_int
IF ( basin_id(sp,sb) < invented_basins ) THEN
DO nb=1,NbNeighb
IF (( neighbours_g(sp,order_ref(nb)) == minhloc(basin_id(sp,sb),1)) .AND. (found == 0)) THEN
dop = minhloc(basin_id(sp,sb),1)
bop = minhloc(basin_id(sp,sb),2)
IF ((dop < undef_int ) .AND. (dop .GT. 0) .AND. (bop < undef_int) .AND. (bop .GT. 0)) THEN
!basin_hierarchy(sp,sb) = basin_hierarchy(dop,bop)
CALL routing_updateflow(sp, sb, dop, bop, nbpt,nwbas, nbxmax*8, outflow_grid, outflow_basin, &
& inflow_number, inflow_grid, inflow_basin)
found = 1
WRITE(numout,*) "Lowest basin hierarchy in the neighbours grid points"
WRITE(numout,*) "hierarch of (sp,sb)", basin_hierarchy(sp,sb)
WRITE(numout,*) "Lowest basinid hierarch", basin_hierarchy(dop,bop)
END IF
END IF
END DO
END IF
ENDIF
!
! Last option : coastal outflow
IF (found .EQ. 0) THEN
!
! Trung: Well, if you come to this step, it should be a special
! (or curious, I said) case. For example, when processing
! Sulina branch of river Danube with the forcing data
! from E2OFD (/bdd/ORCHIDEE_Forcing/BC/OOL/OL2/E2OFD/)
! in resolution of 0.25 degree. The grid box [45.125N, 29.375E]
! is a ocean point in the middle of 7 land points ("a
! hole"). But HydroSHEDS shows that Sulina branch flows
! through this grid cell (from [45.125N, 29.125E] to [45.125N, 29.625E]).
! Trung: If the "JUMP" does not work, it means a worse case.
! Or it means a worse case.
! The outlet of this river falls in the ocean point of
! forcing data. For example: with 0.25 degree forcing
! data of E2OFD, the grid point [ 45.625N ; 13.875E ]
! lost the outlet of basin ID 69666.
! We need to use this point (which should be the last point
! of this river with shortest hierarchy and highest flow accumulation)
! to become new outlet.
!
! This is not good solution because there can be few
! outlet for one basin ID.
! But I need the model works now !!! => so, come back
! here later !
!
WRITE (numout,*) "Linkup : Made a NEW OUTLET at sp & sb: ",sp,sb
! Coastal flow or river flow is both ok here
outflow_grid(sp,sb) = -2
basin_hierarchy(sp,sb) = 0.00
ENDIF
IF ( outflow_basin(sp,sb) < undef_int ) solved(sp,4) = solved(sp,4)+1
ENDIF ENDIF
ENDDO ENDDO
ENDDO ENDDO
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment