SubcatchmentSubarea.java [src/java/m/weppws] Revision: default  Date:
package m.weppws;

import csip.Config;
import csip.api.server.Executable;
import csip.api.client.ModelDataServiceCall;
import csip.api.server.ServiceException;
import csip.SessionLogger;
import csip.utils.Parallel;
import csip.utils.Parallel.Run;
import java.io.BufferedWriter;
import java.io.File;
import java.io.FileWriter;
import java.io.IOException;
import java.util.ArrayList;
import java.util.List;
import java.util.concurrent.atomic.AtomicInteger;
import java.util.function.Function;
import org.codehaus.jettison.json.JSONException;
import org.codehaus.jettison.json.JSONObject;
import util.Grid;
import util.WeppConstants;

/**
 * Represents a WEPP and TauDEM subcatchment. This can be modeled as flowpaths
 * or combined into a representative hillslope.
 *
 */
public class SubcatchmentSubarea {

  int id;
  int weppid;   
  int topazid;
  
  int lastFlowpathID;
  int majorSoil;
  int majorLanduse;
  Subcatchment parSubcatch;
  WeppConstants.SubcatchmentType subType;
  WeppConstants.Orientation dir;
  V1_0 service;
  
  float avgSoilLoss;      // representative soil loss based on all flowpaths
  float avgRunoff;        // representative runoff based on all flowpaths
  float avgSedYield;      // representative sediment yld based on all flowpaths
  
  float repProfileWeights[];
  float repProfileSlopes[];
 
  double area;
  int cells;
  Flowpath repHillslopeFlowpath;    // representative hillslope of all orientations
  double repHillslopeWidth;
  int flowpathsRun;
  int numFlowpaths;
  boolean inUse;
  int weppAltID;


 SubcatchmentSubarea(int id, V1_0 par, int idx, WeppConstants.SubcatchmentType type, Subcatchment par2) {
    this.id = id;
    this.service = par;
    this.parSubcatch = par2;
    this.subType = type;
    weppid = idx;
    repHillslopeFlowpath = null;
    inUse = false;
    majorSoil = par2.majorSoil;    // TBD - this needs to be recalculated
    majorLanduse = par2.majorLanduse;  // TBD - this needs to be recalculated
    area = par2.area; // TBD - this needs to be recalculated
    cells = par2.cells; // TBD - this needs to be recalculated
    weppAltID = 0;
    if (type == WeppConstants.SubcatchmentType.TOP) {
       dir = WeppConstants.Orientation.TOP;
    } else if (type == WeppConstants.SubcatchmentType.LEFT) {
       dir = WeppConstants.Orientation.LEFT;
    } else if (type == WeppConstants.SubcatchmentType.RIGHT) {
       dir = WeppConstants.Orientation.RIGHT;
    }
    topazid = id;
  }
  
 
  int writeFlowpaths(File baseDir) throws IOException, JSONException {
    File fileName = new File(baseDir, "flowpaths_" + Integer.toString(id));
    int num = 0;
    try (BufferedWriter writer = new BufferedWriter(new FileWriter(fileName))) {
      for (int i = 0; i < parSubcatch.flowpaths.size(); i++) {
        if (parSubcatch.flowpaths.get(i).getDir() == dir) {
              parSubcatch.flowpaths.get(i).print(writer);
              num++;
        }
      }
      writer.close();
    }
    return parSubcatch.flowpaths.size();
  }


  String summarizeInput(int hillnum) throws JSONException {
    String ostr = getOrientation();
    
    JSONObject so = new JSONObject()
        .put("id", id)
        .put("weppid", weppid)
        .put("weppAltID", weppAltID)
        .put("topazid", topazid)
        .put("orientation", ostr)
        .put("Area (ac)", area * WeppConstants.CONV_HA_TO_AC)
        .put("Grid cells", cells)
        .put("Major Soil", service.soilIndexToShortName(majorSoil))
        .put("Major Soil id", majorSoil)
        .put("Major Landuse", service.manIndexToShortName(majorLanduse))
        .put("Major Landuse id", majorLanduse)
        .put("Percentage of total area", getPercentageArea());
    return so.toString(2);
  }


  double getPercentageArea() {
      return ((double) cells / (double) service.aoi.getArea()) * 100;
  }
  
  
  String summarizeOutput(SessionLogger log, int hillnum) throws JSONException {
    double avgRunoffVol;
    String ostr = getOrientation();
   
    try {
        avgRunoffVol = (repHillslopeFlowpath.totalLengthMeters*WeppConstants.CONV_M_TO_FT) * 
            (repHillslopeWidth*WeppConstants.CONV_M_TO_FT) * (avgRunoff/12.0);
    } catch (Exception e) {
        avgRunoffVol = -999;
        log.info(e.toString());
        if (repHillslopeFlowpath != null) {
            log.info(repHillslopeFlowpath.toString());
        } else {
            log.info("No representative hillslope created for subcatchment " + id);
            log.info("Number of flowpaths: " + numFlowpaths);
            JSONObject so = new JSONObject()
                .put("id", id)
                .put("weppid", weppid)
                .put("topazid", topazid)
                .put("orientation", ostr)
                .put("Area (ac)", area * WeppConstants.CONV_HA_TO_AC)
                .put("Grid cells", cells)
                .put("Major Soil", service.soilIndexToShortName(majorSoil))
                .put("Major Soil id", majorSoil)
                .put("Major Landuse", service.manIndexToShortName(majorLanduse))
                .put("Major Landuse id", majorLanduse)
                .put("Percentage of total area", ((double) cells / (double) service.aoi.getArea()) * 100)
                .put("Representative Hillslope Length (ft)",-999)
                .put("Representative Hillslope Width (ft)",-999)
                .put("Soil Loss (ton/yr)", -999)
                .put("Soil Loss per Area (ton/ac/yr)", -999)
                .put("Sediment Yield (ton/yr)", -999)
                .put("Sediment Yield per Area (ton/ac/yr)", -999)
                .put("Runoff (in/yr)", -999)
                .put("Runoff Volume (ft^3/yr)", -999)
                .put("Flowpaths within subcatchment", numFlowpaths)
                .put("Message", "No flowpaths calculated in subcatchment");
                return so.toString(2);
        }
    }
    JSONObject so = new JSONObject()
        .put("id", id)
        .put("weppid", weppid)
        .put("topazid", topazid)
        .put("weppAltID", weppAltID)
        .put("orientation", ostr)
        .put("Area (ac)", area * WeppConstants.CONV_HA_TO_AC)
        .put("Grid cells", cells)
        .put("Major Soil", service.soilIndexToShortName(majorSoil))
        .put("Major Soil id", majorSoil)
        .put("Major Landuse", service.manIndexToShortName(majorLanduse))
        .put("Major Landuse id", majorLanduse)
        .put("Percentage of total area", ((double) cells / (double) service.aoi.getArea()) * 100)
        .put("Representative Hillslope Length (ft)",repHillslopeFlowpath.totalLengthMeters*WeppConstants.CONV_M_TO_FT)
        .put("Representative Hillslope Width (ft)",repHillslopeWidth*WeppConstants.CONV_M_TO_FT)
        .put("Soil Loss (ton/yr)", avgSoilLoss)
        .put("Soil Loss per Area (ton/ac/yr)", avgSoilLoss/(area * WeppConstants.CONV_HA_TO_AC))
        .put("Sediment Yield (ton/yr)", avgSedYield)
        .put("Sediment Yield per Area (ton/ac/yr)", avgSedYield/(area * WeppConstants.CONV_HA_TO_AC))
        .put("Runoff (in/yr)", avgRunoff)
        .put("Runoff Volume (ft^3/yr)", avgRunoffVol)
        .put("Flowpaths within subcatchment", numFlowpaths)
        .put("Message","");
    return so.toString(2);
  }


  //
  // runFlowpaths0()
  //
  // Uses CSIP services to make parallel WEPP runs.
  //
  int runFlowpaths0(SessionLogger LOG, File workspace, Grid field, boolean statsOnly) throws Exception {

    LOG.info("Parallel Flowpath runs for subcatchment " + id);
    AtomicInteger numRun = new AtomicInteger();
    int nthreads = Config.getInt("weppws.fp.par", 16);

    Parallel.run(nthreads, parSubcatch.flowpaths.size(), new Function<Integer, Run>() {
      @Override
      public Run apply(Integer t) {
        System.out.println("Created FP Run: " + t);
        service.progress("Running " + id + "/" + t);
        return () -> {
          Flowpath fp = parSubcatch.flowpaths.get(t);
          if (!statsOnly && fp.intersectsField(field) && parSubcatch.flowpaths.get(t).getDir() == dir) {
            System.out.println("Running FP: >>>>>>>>>>>>>>>>>  " + t);
            fp.runFlowpathService(workspace, service.aoi, service.flow_d8_grid.getCellsize(),false);
            numRun.incrementAndGet();
          }
        };
      }
    });

    return flowpathsRun = numRun.get();
  }

  //
  // runFlowpaths1()
  //
  // Runs all, or a subset percentage, of the flowpaths in a the subcatchment.
  // The runs are done by calling the local executable, a service call is not 
  // done to run WEPP. Every call must wait until the previous run is complete.
  //
  int runFlowpaths1(SessionLogger LOG, File workspace, Executable weppserv, Grid field,
      boolean statsOnly, int subsetPercentage) throws ServiceException, IOException, Exception {
      
    LOG.info("Exe Flowpath runs for subcatchment " + id);
    int numRun = 0;
    double div = ((double) subsetPercentage) / (double) 100;
    int freq = (int) (1.0 / div);

    flowpathsRun = 0;
    int i = 0;
    while (i < parSubcatch.flowpaths.size()) {
      if (parSubcatch.flowpaths.get(i).intersectsField(field) && parSubcatch.flowpaths.get(i).getDir() == dir) {
        LOG.info(">>> Starting flowpath: " + i);
        if (statsOnly == false) {
          // may only be running a subset of the flowpaths
          if ((subsetPercentage == 100) || ((i % freq) == 0)) {
              parSubcatch.flowpaths.get(i).runFlowpathLocal(workspace, weppserv,
                  service.aoi, service.flow_d8_grid.getCellsize(), false);
              numRun++;
          } 
        } else {
          numRun++;
        }
      } 
      i++;
    }
    flowpathsRun = numRun;
    return numRun;
  }
  
  

  void makeRepresentativeHillslope(float cellSize) {
    float longest = 0;
    int matched = 0;
    for (int i = 0; i < parSubcatch.flowpaths.size(); i++) {
      if (parSubcatch.flowpaths.get(i).getDir() == dir) {
          matched++;
          if (parSubcatch.flowpaths.get(i).totalLengthMeters > longest) {
            longest = parSubcatch.flowpaths.get(i).totalLengthMeters;
          }
      }
    }
    //log.info("  ... Matched: " + matched + " out of: " + parSubcatch.flowpaths.size() + " Longest: " + longest);
    int cellsNeeded = (int) ((double) longest / cellSize);
    repProfileWeights = new float[cellsNeeded];
    repProfileSlopes = new float[cellsNeeded];

    numFlowpaths = 0;
    
    for (int i = 0; i < parSubcatch.flowpaths.size(); i++) {
      if (parSubcatch.flowpaths.get(i).getDir() == dir) {
         parSubcatch.flowpaths.get(i).updateRepProfile(cellSize, repProfileSlopes, repProfileWeights);
         numFlowpaths++;
      }
    }
    for (int i = 0; i < repProfileSlopes.length; i++) {
      repProfileSlopes[i] = repProfileSlopes[i] / repProfileWeights[i];
    }

    int r[] = new int[cellsNeeded];
    int c[] = new int[cellsNeeded];
    float sl[] = new float[cellsNeeded];
    int s[] = new int[cellsNeeded];
    int m[] = new int[cellsNeeded];
    float d[] = new float[cellsNeeded];

    for (int i = 0; i < cellsNeeded; i++) {
      r[i] = 0;
      c[i] = 0;
      sl[i] = repProfileSlopes[i] / (float) 100.0;   // slopes were calculated as %, need to be 0-1 for constructor
      s[i] = majorSoil;
      m[i] = majorLanduse;
      d[i] = (float) cellSize;
    }
    repHillslopeFlowpath = new Flowpath(topazid, r, c, sl, m, s, d, service, -1, -1, -1);
    repHillslopeWidth = (area * 10000) / repHillslopeFlowpath.totalLengthMeters;
  }


 
        
  int runRepresentativeHillslope0(File workspace, Executable weppserv,
      Grid field, boolean statsOnly, boolean usePassFiles) throws ServiceException, IOException, Exception {
      
    makeRepresentativeHillslope((float) field.getCellsize());
    
    if (repHillslopeFlowpath == null) {
        avgSoilLoss = 0;
        avgRunoff = 0;
        avgSedYield = 0;
        return 0;
    }
    
    if (statsOnly == false) {
      repHillslopeFlowpath.runFlowpathService(workspace, service.aoi, repHillslopeWidth, usePassFiles);
      avgSoilLoss = repHillslopeFlowpath.avgSoilLoss;
      avgRunoff = repHillslopeFlowpath.avgRunoff;
      avgSedYield = repHillslopeFlowpath.avgSedYield;
    }
    return 1;
  }
  
  // 
  // runRepresentativeHillslope()
  //
  // Runs the flowpath that is an aggregation of all flowpaths in a subcatchment.
  // Calls to either run the simulation as a service or directly (runflowpath1).
  //
   JSONObject runRepresentativeHillslope(SessionLogger log, File workspace, Executable weppserv,
      Grid field, boolean statsOnly, boolean usePassFiles) throws ServiceException, IOException, Exception {
      
    JSONObject resp = null;
   
    makeRepresentativeHillslope((float) field.getCellsize());
    
    //log.info("   ...runRep Length is: " + repHillslopeFlowpath.totalLengthMeters);
    
    if (repHillslopeFlowpath == null) {
        avgSoilLoss = 0;
        avgRunoff = 0;
        avgSedYield = 0;
        return null;
    }
    if (statsOnly == false) {
      //if (service.flowpathURL.length() > 0) {
      //     resp = repHillslopeFlowpath.runFlowpath(workspace, weppserv, service.aoi, repHillslopeWidth);
      //     if (resp != null) {
      //       avgSoilLoss = repHillslopeFlowpath.avgSoilLoss;
      //        avgRunoff = repHillslopeFlowpath.avgRunoff;
      //        avgSedYield = repHillslopeFlowpath.avgSedYield;
      //      } 
      //} else {
           resp = repHillslopeFlowpath.runFlowpathLocal(workspace, weppserv, service.aoi, repHillslopeWidth, usePassFiles);
           avgSoilLoss = repHillslopeFlowpath.avgSoilLoss;
           avgRunoff = repHillslopeFlowpath.avgRunoff;
           avgSedYield = repHillslopeFlowpath.avgSedYield;
      //}
     
    }
    return resp;
  }

 //
 // checkProgress()
 //
 // Checksif the representaive hillslope run is complete if it was started from a
 // service call.
 //
   JSONObject checkProgress(SessionLogger LOG, JSONObject resp, File workspace) 
        throws ServiceException, IOException, Exception {
       JSONObject resp2 = repHillslopeFlowpath.checkProgress(LOG, resp, workspace, service.aoi);
       if (resp2 == null) {
            avgSoilLoss = repHillslopeFlowpath.avgSoilLoss;
            avgRunoff = repHillslopeFlowpath.avgRunoff;
            avgSedYield = repHillslopeFlowpath.avgSedYield;
       }
       return resp2;
   }
           
  float getAvgSoilLoss() {
    return avgSoilLoss;
  }


  float getAvgRunoff() {
    return avgRunoff;
  }


  float getAvgSedYield() {
    return avgSedYield;
  }


  Grid outputFlowpaths(Grid lossGrid, Grid sedWeightGrid) {
    for (int i = 0; i < parSubcatch.flowpaths.size(); i++) {
      if (parSubcatch.flowpaths.get(i).getDir() == dir) {  
         parSubcatch.flowpaths.get(i).updateWeightedLossGrid(lossGrid, sedWeightGrid);
      }
    }
    return lossGrid;
  }


  void setMajorSoil(int val) {
    majorSoil = val;
  }


  void setMajorLanduse(int val) {
    majorLanduse = val;
  }


  void calcArea(Grid subcatchments, Grid bound) throws ServiceException {
    cells = subcatchments.cellCountWithMask(topazid, bound);
    area = (cells * subcatchments.getCellsize() * subcatchments.getCellsize()) / 10000; // total area in ha
  }

  String printProfile() {
    if (repHillslopeFlowpath != null) {
        return repHillslopeFlowpath.printProfile(repHillslopeWidth);
    } else {
        return("Representative profile not set for hillslope." + String.valueOf(id));
    }
  }

String soilIndexToName(int idx) throws ServiceException {
    String name;
    try {
      name = service.soilIndexToName(idx);
      if (name.equals("")) {
        throw new ServiceException("Could not find soil to match index: " + Integer.toString(idx));
      }
    } catch (JSONException je) {
      throw new ServiceException("Could not find soil to match index: " + Integer.toString(idx));
    }
    return name;
  }


  String manIndexToName(int idx) throws ServiceException {
    String name;
    try {
      name = service.manIndexToName(idx);
      if (name.equals("")) {
        throw new ServiceException("Could not find management to match index: " + Integer.toString(idx));
      }
    } catch (JSONException je) {
      throw new ServiceException("Could not find management to match index: " + Integer.toString(idx));
    }
    return name;
  }

  double getRepHillslopeWidth() {
    return repHillslopeWidth;
  }


  double getRepHillslopeLength() {
    return repHillslopeFlowpath.totalLengthMeters;
  }


  String majorLanduseName() throws ServiceException {
    return manIndexToName(majorLanduse);
  }


  String majorSoilName() throws ServiceException {
    return soilIndexToName(majorSoil);
  }
  
  void setInuse(boolean val) {
      inUse = val;
  }
 
  boolean getInuse() {
      return inUse;
  }
  
  void setWeppID(int val) {
      weppid = val;
  }
  
   void setWeppAltID(int val) {
      weppAltID = val;
  }
   
  int getWeppAltID() {
      return weppAltID;
  }
  int getCells() {
      return cells;
  }
  
  String getOrientation() {
    String ostr;
    if (dir == WeppConstants.Orientation.LEFT) {
        ostr = "LEFT";
    } else if (dir == WeppConstants.Orientation.RIGHT) {
        ostr = "RIGHT";
    } else if (dir == WeppConstants.Orientation.TOP) {
        ostr = "TOP";
    } else {
        ostr = "UNKNOWN";
    }
    return ostr;
  }
}