|
SURFEX v7.3
General documentation of Surfex
|
|
SURFEX v7.3
General documentation of Surfex
|
00001 !----------------------------------------------------------------- 00002 ! ####################### 00003 SUBROUTINE MAKE_MASK_TOPD_TO_ISBA(KI) 00004 ! ####################### 00005 ! 00006 !!**** *MAKE_MASK_TOPD_TO_ISBA(* 00007 !! 00008 !! PURPOSE 00009 !! ------- 00010 ! 00011 ! Create a mask for each catchment. 00012 ! 00013 ! 00014 !!** METHOD 00015 !! ------ 00016 ! 00017 !! EXTERNAL 00018 !! -------- 00019 !! 00020 !! none 00021 !! 00022 !! IMPLICIT ARGUMENTS 00023 !! ------------------ 00024 !! 00025 !! REFERENCE 00026 !! --------- 00027 00028 !! AUTHOR 00029 !! ------ 00030 !! 00031 !! L. Labatut & K. Chancibault * CNRM * 00032 !! 00033 !! MODIFICATIONS 00034 !! ------------- 00035 !! 00036 !! Original 16/03/2005 00037 !------------------------------------------------------------------------------- 00038 ! 00039 !* 0. DECLARATIONS 00040 ! ------------ 00041 ! 00042 USE MODD_TOPODYN, ONLY: NNCAT, NNYC, XY0, XDXT, NNXC,& 00043 XX0, XTOPD, XNUL, NLINE, NMESHT 00044 USE MODD_COUPLING_TOPD, ONLY: NMASKT 00045 USE MODD_SURF_PAR, ONLY : XUNDEF, NUNDEF 00046 ! 00047 USE MODI_GET_LUOUT 00048 USE MODI_ABOR1_SFX 00049 USE MODI_WRITE_FILE_MAP 00050 ! 00051 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00052 USE PARKIND1 ,ONLY : JPRB 00053 ! 00054 IMPLICIT NONE 00055 ! 00056 !* 0.1 declarations of arguments 00057 ! 00058 INTEGER, INTENT(IN) :: KI ! Grid dimensions 00059 ! 00060 !* 0.2 declarations of local variables 00061 ! 00062 CHARACTER(LEN=30) :: YVAR ! name of results file 00063 INTEGER :: JCAT, JJ, JI, IDX ! loop control 00064 INTEGER :: II,ILINE ! work integer variables 00065 INTEGER :: IDXM ! indexes of Isba grid meshes and nodes 00066 INTEGER :: ILUOUT ! unit 00067 REAL :: ZXT, ZYT ! catchment grid nodes Lambert II coordinates 00068 REAL :: ZX1, ZX2, ZX3, ZX4, ZY1, ZY2, ZY3, ZY4 ! Isba mesh Lambert II coordinates 00069 REAL :: ZXA, ZXB, ZYA, ZYB 00070 REAL, DIMENSION(NNCAT,NMESHT):: ZWRK 00071 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00072 !------------------------------------------------------------------------------- 00073 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA',0,ZHOOK_HANDLE) 00074 ! 00075 CALL GET_LUOUT('OFFLIN',ILUOUT) 00076 ! 00077 !WRITE(*,*) 'Y0=', XY0(1), XY0(2) 00078 !WRITE(*,*) 'X0=', XX0(1), XX0(2) 00079 !WRITE(*,*) 'X et Y isba=', XXI(1), XXI(180), XYI(1), XYI(180) 00080 !loop on catchments 00081 DO JCAT=1,NNCAT 00082 ! 00083 IDXM = 1 00084 CALL INIT_4POINTS(IDXM,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4) 00085 ! 00086 DO JJ=1,NNYC(JCAT) ! number of topographic points on the y axis 00087 ! 00088 ZYT = XY0(JCAT) + (JJ-1) * XDXT(JCAT) 00089 ZYT = ZYT + 0.5 * XDXT(JCAT) 00090 ! 00091 DO JI=1,NNXC(JCAT) 00092 ! 00093 ZXT = XX0(JCAT) + (JI-1) * XDXT(JCAT) 00094 ZXT = ZXT + 0.5 * XDXT(JCAT) 00095 ! 00096 IDX = (JJ-1) * NNXC(JCAT) + JI ! index of the point among all in the catchment 00097 ! 00098 !* on vérifie que le pixel du MNT appartient au BV 00099 IF ( XTOPD(JCAT,IDX).NE.XNUL(JCAT) ) THEN 00100 !* calcule des coordonées X et Y inf et sup de la maille ISBA considérée 00101 CALL GET_COORD(ZXT,ZYT,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4,ZXA,ZYA,ZXB,ZYB) 00102 !* si on se trouve sur le premier pixel du MNT ou si le pixel du MNT n'est pas 00103 !dans la maille Isba considérée (qui est celle dans laquelle se trouve le pixel précédent) 00104 IF (ZXT.LT.ZXA.OR.ZXT.GE.ZXB.OR.ZYT.LT.ZYA.OR.ZYT.GE.ZYB) THEN 00105 !* on repart de la première maille de la grille Isba 00106 IDXM = 1 00107 CALL INIT_4POINTS(IDXM,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4) 00108 CALL GET_COORD(ZXT,ZYT,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4,ZXA,ZYA,ZXB,ZYB) 00109 !* on parcours les mailles de la grille Isba, jusqu'à ce qu'on trouve la maille à laquelle appartient le pixel du MNT 00110 DO WHILE (ZXT.LT.ZXA.OR.ZXT.GE.ZXB.OR.ZYT.LT.ZYA.OR.ZYT.GE.ZYB) 00111 IDXM = IDXM + 1 00112 IF (IDXM.GE.KI) THEN 00113 WRITE(*,*) 'ZXT', ZXT,'ZYT',ZYT 00114 WRITE(*,*) 'indices Isba:',IDXM,'>=',KI 00115 CALL ABOR1_SFX("MAKE_MASK_TOPD_TO_ISBA: PROBLEM") 00116 ENDIF 00117 CALL INIT_4POINTS(IDXM,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4) 00118 CALL GET_COORD(ZXT,ZYT,ZX1,ZX2,ZX3,ZX4,ZY1,ZY2,ZY3,ZY4,ZXA,ZYA,ZXB,ZYB) 00119 ENDDO 00120 ENDIF 00121 IF (NLINE(JCAT,IDX)/=0) NMASKT(JCAT,NLINE(JCAT,IDX)) = IDXM 00122 ENDIF 00123 ENDDO 00124 ENDDO 00125 ENDDO 00126 ! 00127 YVAR='.mask_topd' 00128 WHERE (NMASKT(:,:)/=NUNDEF) 00129 ZWRK(:,:)=REAL(NMASKT) 00130 ELSEWHERE 00131 ZWRK(:,:)=XUNDEF 00132 ENDWHERE 00133 CALL WRITE_FILE_MAP(ZWRK,YVAR) 00134 !CALL WRITE_FILE_MAP(REAL(NMASKT),YVAR) 00135 ! 00136 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA',1,ZHOOK_HANDLE) 00137 ! 00138 CONTAINS 00139 ! 00140 SUBROUTINE INIT_4POINTS(KDXM,PX1,PX2,PX3,PX4,PY1,PY2,PY3,PY4) 00141 ! 00142 USE MODD_COUPLING_TOPD, ONLY: NIMAX, XXI, XYI 00143 ! 00144 INTEGER, INTENT(IN) :: KDXM 00145 REAL, INTENT(OUT) :: PX1, PX2, PX3, PX4 00146 REAL, INTENT(OUT) :: PY1, PY2, PY3, PY4 00147 ! 00148 INTEGER :: ILINE, II, IDXN 00149 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00150 ! 00151 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:INIT_4POINTS',0,ZHOOK_HANDLE) 00152 ! 00153 ILINE = INT(KDXM/(NIMAX))+1 ! number of the current line 00154 II = KDXM-((ILINE-1)*NIMAX) ! index of point in the line 00155 IDXN = (ILINE-1)*(NIMAX+1)+II ! indice du point dans la grille isba 00156 ! 00157 PX1 = XXI(IDXN) ! coordonnée X du point courant 00158 PX2 = XXI(IDXN+1) ! coordonnée X du point suivant 00159 PX3 = XXI(IDXN+(NIMAX+1)) ! coordonnée X du point aligné sur la ligne suivante 00160 PX4 = XXI(IDXN+1+(NIMAX+1)) ! coordonnée X du point suivant sur la ligne suivante 00161 ! 00162 PY1 = XYI(IDXN) ! coordonnée Y du point courant 00163 PY2 = XYI(IDXN+1) ! coordonnée Y du point suivant 00164 PY3 = XYI(IDXN+(NIMAX+1)) ! coordonnée Y du point aligné sur la ligne suivante 00165 PY4 = XYI(IDXN+1+(NIMAX+1)) ! coordonnée Y du point suivant sur la ligne suivante 00166 ! 00167 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:INIT_4POINTS',1,ZHOOK_HANDLE) 00168 ! 00169 END SUBROUTINE INIT_4POINTS 00170 ! 00171 SUBROUTINE GET_COORD(PXT,PYT,PX1,PX2,PX3,PX4,PY1,PY2,PY3,PY4,& 00172 PXA,PYA,PXB,PYB) 00173 ! 00174 REAL, INTENT(IN) :: PXT, PYT 00175 REAL, INTENT(IN) :: PX1, PX2, PX3, PX4 00176 REAL, INTENT(IN) :: PY1, PY2, PY3, PY4 00177 REAL, INTENT(OUT) :: PXA, PYA, PXB, PYB 00178 REAL :: ZFA, ZFB 00179 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00180 ! 00181 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:GET_COORD',0,ZHOOK_HANDLE) 00182 ! 00183 IF((PX3-PX1).EQ.0.0) THEN 00184 PXA = PX1 00185 ELSE 00186 CALL GET_LINE_PARAM(PX1,PY1,PX3,PY3,ZFA,ZFB) 00187 PXA = (ZYT-ZFB)/ZFA 00188 ENDIF 00189 ! 00190 IF ((PX4-PX2).EQ.0.0) THEN 00191 PXB = PX2 00192 ELSE 00193 CALL GET_LINE_PARAM(PX2,PY2,PX4,PY4,ZFA,ZFB) 00194 PXB = (ZYT-ZFB)/ZFA 00195 ENDIF 00196 ! 00197 IF ((PY2-PY1).EQ.0.0) THEN 00198 PYA = PY2 00199 ELSE 00200 CALL GET_LINE_PARAM(PX1,PY1,PX2,PY2,ZFA,ZFB) 00201 PYA = ZFA*ZXT+ZFB 00202 ENDIF 00203 ! 00204 IF ((PY4-PY3).EQ.0.0) THEN 00205 PYB = PY4 00206 ELSE 00207 CALL GET_LINE_PARAM(PX3,PY3,PX4,PY4,ZFA,ZFB) 00208 PYB = ZFA*ZXT+ZFB 00209 ENDIF 00210 ! 00211 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:GET_COORD',1,ZHOOK_HANDLE) 00212 ! 00213 END SUBROUTINE GET_COORD 00214 ! 00215 SUBROUTINE GET_LINE_PARAM(PX1,PY1,PX2,PY2,PFA,PFB) 00216 ! 00217 REAL, INTENT(IN) :: PX1, PX2, PY1, PY2 00218 REAL, INTENT (OUT) :: PFA, PFB 00219 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00220 ! 00221 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:GET_LINE_PARAM',0,ZHOOK_HANDLE) 00222 ! 00223 PFA = (PY2 - PY1) / (PX2 - PX1) ! slope 00224 PFB = PY1 - PFA * PX1 ! offset 00225 ! 00226 IF (LHOOK) CALL DR_HOOK('MAKE_MASK_TOPD_TO_ISBA:GET_LINE_PARAM',1,ZHOOK_HANDLE) 00227 ! 00228 END SUBROUTINE GET_LINE_PARAM 00229 ! 00230 END SUBROUTINE MAKE_MASK_TOPD_TO_ISBA
1.8.0
