      PROGRAM MAIN_PROG1

      IMPLICIT NONE

      ! THIS IS THE MAIN PROGRAM FROM WHERE ALL THE OTHER SUBROUTINES
      ! ARE CALLED. THIS PROGRAM ACCEPTS THE DAILY MELTWATER RUNOFF,
      ! DAILY SNOW DEPTH, DAILY SNOW DENSITY AND THE DEM OF THE STUDY
      ! AREA

      INTEGER, PARAMETER :: RO=486,CO=3104,SX=10000

      INTEGER :: ROWS,COLS,FLINDI(CO,RO),FLINDJ(CO,RO),K,
     &CMARK(CO,RO),FLAG_SINK(CO,RO),OV_LAKE(CO,RO),
     &SINK_I(SX),SINK_J(SX),LAKE_I(SX),LAKE_J(SX),
     &MARK2(CO,RO),DRAIN_INCISE(CO,RO),LAK_NO_INC(CO,RO),
     &DRAIN(CO,RO),NUC_I(SX),NUC_J(SX),SSX,I,J,FLAG_OV(CO,RO)

      DOUBLE PRECISION :: SURF(CO,RO),DEM(CO,RO),
     &RUNOFF(CO,RO),NEW_RU(CO,RO),SDEP(CO,RO),
     &SDEN(CO,RO),STRESS(CO,RO),ICE(CO,RO),
     &CREV_DEP(CO,RO),WAT_DEP(CO,RO),LAKE_DEPTH(CO,RO),
     &LAKE_HT(CO,RO),LAKE_AR(CO,RO),LAKE_BASE(CO,RO),
     &LAKE_RED(CO,RO),CHN_DEC(CO,RO),CHN_DEC1(CO,RO),
     &RUN_X(CO,RO),RUN_Y(CO,RO),RUN_LAK(CO,RO),
     &LAKE2RUN(CO,RO),TOTAL_DRAIN(CO,RO),SUMINC(CO,RO),
     &MELT(CO,RO),RAIN(CO,RO),SUML(CO,RO)
      

      CHARACTER :: PATH*42,FILENAME*53,FILEN*54,DCHAR*4,FNAME*60,
     &FNAME1*60,FNAME2*60,FNAME3*56

      EXTERNAL HYDROFRACTURE1,
     &CREV_LAKE1,
     &STRESS_INT1,
     &NON_LAKE_CREV_NEW1,
     &CREV_CELLS1,
     &LAKE_CHN_INC_OV1, 
     &LAKE_DEL1,
     &RUNOFF_DIST1,
     &LAKE_CATCHMENT1,
     &NEW_OUTLET1, 
     &PATH_SEL_CHN_INC1,
     &CHANNEL_INC1,
     &ONLY_OVERFLOW1,
     &LAKE_RED_TRANS1,
     &FLOW_DIR1,
     &WATER_TRANSFER1,
     &SINK_SEPARATE1,
     &FILL_DEPRESSION1,
     &FILL_HALF_FULL_DEPRESSION1,
     &MARGIN_FLOW1,
     &LAKE_NUCLEUS_ISOLATION1,
     &SINK_LOCATOR1,
     &LAKE_NUCLEUS_ISOLATION21,
     &REINIT1 

      ! SPECIFYING DIMENSIONS OF THE STUDY AREA I.E. NUMBER OF ROWS AND
      ! COLUMNS

      ROWS=RO
      COLS=CO
      SSX=SX

      ! INITIALISING THE VARIABLES

      FLINDI=0.
      FLINDJ=0.
      SUML=0.
      ICE=0.
      CREV_DEP=0.
      WAT_DEP=0.
      FLAG_SINK=0.
      LAKE_DEPTH=0.
      OV_LAKE=0.
      LAKE_HT=0.
      LAKE_AR=0.
      LAKE_BASE=0.
      LAKE_RED=0.
      CHN_DEC=0.
      SUMINC=0.
      CHN_DEC1=0.
      RUN_X=0.
      RUN_Y=0.
      RUN_LAK=0.
      DRAIN_INCISE=0.
      LAKE2RUN=0.
      LAK_NO_INC=0.
      SINK_I=0.
      SINK_J=0.
      MELT=0.
      RAIN=0.
      STRESS=0.
      CMARK=0.

      PATH='/mnt2/prateek/MAR_data_paper/paper_basin7/'

      ! READ DEM

      OPEN(UNIT=1,FILE=PATH//'dem2.dat',STATUS='OLD')
      READ(1,*) DEM
     
      SURF=DEM

      ! READ ICE THICKNESS

      OPEN(UNIT=15,FILE=PATH//'Ice_thickness.dat',STATUS='OLD')
      READ(15,*) ICE

      DO 100 K=100,334

      ! READ DAILY RUNOFF

      WRITE (DCHAR, '(I4.4)') K
      FILENAME=PATH//'day'//DCHAR//'.dat'
      OPEN(UNIT=1,FILE=FILENAME,STATUS='OLD')
      READ(1,*) NEW_RU
      NEW_RU=NEW_RU/1000. ! IN M/DAY

      ! READ DAILY SNOW DEPTH

      FILEN=PATH//'sdep'//DCHAR//'.dat'
      OPEN(UNIT=11,FILE=FILEN,STATUS='OLD')
      READ(11,*) SDEP

      ! READ DAILY SNOW DENSITY

      FILEN=PATH//'sden'//DCHAR//'.dat'
      OPEN(UNIT=12,FILE=FILEN,STATUS='OLD')
      READ(12,*) SDEN

      ! READ DAILY STRESS

      FNAME3=PATH//'stress'//DCHAR//'.dat'
      OPEN(UNIT=21,FILE=FNAME3,STATUS='OLD')
      READ(21,*) STRESS
      STRESS=STRESS/1000.

      ! IDENTIFYING THE CREVASSED CELLS AND SETTING THE CREVASSE DEPTH

      CALL CREV_CELLS1(ROWS,COLS,CMARK,CREV_DEP,STRESS,ICE)


      CHN_DEC1=CHN_DEC1+CHN_DEC ! TO FIND THE TOTAL DEPTH OF CHANNEL INCISED
      OV_LAKE=0. ! CLEAR THE OVERFLOW LAKE MATRIX
      RUNOFF=RUNOFF+NEW_RU ! UPDATE RUNOFF WITH DAILY PRODUCED MELT

      ! CLEARING THE MATRICES THAT WERE USED IN SIMULATING
      ! CHANNEL INCISION AND HYDROFRACTURE

      RUN_X=0.
      RUN_Y=0.
      RUN_LAK=0.
      SUMINC=0.
      LAKE_RED=0.
      LAK_NO_INC=0.
      TOTAL_DRAIN=0.
      FLAG_SINK=0.

      ! CLEARING THE 1D ARRAYS

      SINK_I=0.
      SINK_J=0.
      NUC_I=0.
      NUC_J=0.

      ! STORING THOSE LAKES THAT DRAINED RAPIDLY OR OVERTOPPED

      FNAME2=PATH//'Wate_depth'//DCHAR//'.dat'
      OPEN(UNIT=22,FILE=FNAME2,STATUS='NEW')
      WRITE(22,*) DRAIN

      FNAME2=PATH//'Crev_depth'//DCHAR//'.dat'
      OPEN(UNIT=23,FILE=FNAME2,STATUS='NEW')
      WRITE(23,*) FLAG_OV

      CHN_DEC=0.

      ! IF THERE IS WATER WHICH IS FREE TO FLOW, RUN THE MODEL

      IF (SUM(RUNOFF).GT.0) THEN

              ! LOCATING THE POTENTIAL SINKS IN THE DEM

              CALL SINK_LOCATOR1(SURF,DEM,COLS,ROWS,FLAG_SINK)

              ! FINDING DESTINATION CELL FOR EVERY NON-LAKE CELL

              CALL FLOW_DIR1(SURF,FLINDI,FLINDJ,ROWS,COLS,CMARK,SDEP,
     &SDEN,DEM,RUNOFF,SSX,ICE)

              ! TRANSFERRING WATER FROM SOURCE TO DESTINATIO CELLS

              CALL WATER_TRANSFER1(ROWS,COLS,FLINDI,FLINDJ,RUNOFF)

              ! SEPARATING THE SINKS FROM REST OF DEM PIXELS

              CALL SINK_SEPARATE1(FLAG_SINK,ROWS,COLS,SINK_I,SINK_J,
     &SURF,DEM,SSX)

              ! FILLING EMPTY DEPRESSIONS

              CALL FILL_DEPRESSION1(ROWS,COLS,DEM,RUNOFF,CMARK
     &,SURF,SSX,SINK_I,SINK_J)

              SINK_I=0.
              SINK_J=0.

              ! SEPARATING LAKE NUCLEUS OF FILLED LAKES

              CALL LAKE_NUCLEUS_ISOLATION1(SURF,DEM,NUC_I,NUC_J,ROWS,
     &COLS,SSX,RUNOFF)

              ! FILLING NON-EMPTY DEPRESSIONS

              IF (COUNT(NUC_I.NE.0).GT.0) THEN

                      CALL FILL_HALF_FULL_DEPRESSION1(ROWS,COLS,SURF,
     &DEM,RUNOFF,CMARK,NUC_I,NUC_J,SSX)

              END IF

              NUC_I=0.
              NUC_J=0.

              ! SEPARATING THE LAKE NUCLEUS FROM REST OF THE DEM CELLS

              CALL LAKE_NUCLEUS_ISOLATION21(SURF,DEM,NUC_I,NUC_J,ROWS
     &,COLS,SSX)

              ! SIMULATING NON-LAKE AND  LAKE HYDROFRACTURE

              CALL HYDROFRACTURE1(CMARK,DRAIN,ICE,STRESS,ROWS,COLS,
     &WAT_DEP,SURF,DEM,RUNOFF,CREV_DEP,TOTAL_DRAIN,SSX,NUC_I,NUC_J,SUML)

              NUC_I=0.
              NUC_J=0.

              ! SEPARATING THE LAKE NUCLEUS FROM REST OF THE DEM CELLS

               CALL LAKE_NUCLEUS_ISOLATION21(SURF,DEM,NUC_I,NUC_J,ROWS
     &,COLS,SSX)

              ! SIMULATING CHANNEL INCISION 

              IF (COUNT(NUC_I.NE.0).GT.0) THEN

              ! CHANNEL INCISION ROUTINE

                      CALL LAKE_CHN_INC_OV1(NUC_I,NUC_J,SURF,DEM,ICE,
     &CMARK,SUMINC,LAKE_HT,LAKE_AR,LAKE_BASE,CHN_DEC,RUN_X,RUN_Y,
     &RUN_LAK,DRAIN_INCISE,LAKE_RED,RUNOFF,SDEP,SDEN,LAKE2RUN,
     &LAK_NO_INC,ROWS,COLS,SSX,FLAG_OV)

             ! UPDATING LAKE AND RUNOFF DEPTH AFTER CHANNEL INCISION

                      CALL LAKE_RED_TRANS1(RUNOFF,RUN_LAK,RUN_X,RUN_Y,
     &SURF,DEM,CHN_DEC,LAKE_RED,LAKE2RUN,COLS,ROWS,LAK_NO_INC)

              END IF

              NUC_I=0.
              NUC_J=0.

      END IF

      ! CLEAR RUNOFF FROM THE MARGIN CELLS

      CALL MARGIN_FLOW1(RUNOFF,COLS,ROWS,ICE,DEM)

      LAKE_DEPTH=SURF-DEM

      ! STORING LAKE DEPTH FOR THE DAY

      FNAME=PATH//'Lake_depth'//DCHAR//'.dat'
      OPEN(UNIT=13, FILE=FNAME, STATUS='NEW')
      WRITE(13,*) LAKE_DEPTH

      ! STORING LAKES THAT DRAINED VIA  CHANNEL INCISION

      FNAME1=PATH//'Runo_depth'//DCHAR//'.dat'
      OPEN(UNIT=20,FILE=FNAME1,STATUS='NEW')
      WRITE(20,*) DRAIN_INCISE

      ! STORING DEPTHS OF INCISED CHANNELS

      FNAME1=PATH//'Chan_depth'//DCHAR//'.dat'
      OPEN(UNIT=20,FILE=FNAME1,STATUS='NEW')
      WRITE(20,*) CHN_DEC1

      PRINT *,K,SUM(TOTAL_DRAIN),SUM(SUMINC),SUM(SUML)
      SUML=0.

      ! REINITIALISE THE RUNOFF AND CREVASSED CELL FLAG I.E. CMARK FOR
      ! THE NEXT YEAR

      IF ((K.EQ.365).OR.(K.EQ.730).OR.(K.EQ.1095).OR.(K.EQ.1461)
     &.OR.(K.EQ.1826).OR.(K.EQ.2191).OR.(K.EQ.2556))THEN
                 CALL REINIT1(RO,CO,CREV_DEP,WAT_DEP,ICE,RUNOFF,
     &CMARK)
                 DRAIN=0.
                 DRAIN_INCISE=0.
      END IF

100   CONTINUE

      END 
