#coding=windows-1252
#coding=windows-1252
# This first line, #coding=windows-1252, is required to run the program.

# Purpose: to graph in Python the surface contour plot of average air temperature metric across 11 values of tree cover TC and impervious cover IC, a total of 121 unique IC and TC combinations
# Requires: i-Tree Cool Air model output in the form of Blockgroup results (daily and hourly), with a unique block group for each of the 121 unique TC and IC combinations.
# Requires: Seems to need at least 5 days of outputs to construct averages
# At command prompt, give the path to the folder containing the Blockgroup results
# Example of the command prompt C:\iTree\projects>python z_Map_CoolAirTair_11x11_TC-IC_combinations.py C:\iTree\UnifiedHydro\TestingFilesAndScript\TestCases\SpatialCoolAir\LCScale_11x11Grid\expectedoutput

# Created by ESF student Li Zhang, working with Professor T. Endreny in spring 2020
# V4 11/18/2020

# Requires python 3 interpreter to use matplotlib (and associated mpl_toolkits)
from mpl_toolkits import mplot3d
from mpl_toolkits.mplot3d import Axes3D
import array
import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter
import numpy as np
import linecache
import pandas as pd
import os
from os import path
import sys
# from matplotlib.mlab import griddata #has been removed
from scipy.interpolate import griddata

mod = sys.modules[__name__]

#file path to model inputs, read from command line >python z_Map_CoolAirTair_11x11_TC-IC_combinations.py C:\iTree\UnifiedHydro\TestingFilesAndScript\TestCases\SpatialCoolAir\LCScale_11x11Grid\expectedoutput
#backslash is added to path input
fpath = str(sys.argv[1]  + '\\')
#sysag = "C:\iTree\HydroPlus\TestingFilesAndScript\TestCases\SpatialCoolAir\LCScale_11x11Grid\output"
#sysag = "C:\iTree\HydroPlus\TestingFilesAndScript\TestCases\SpatialCoolAir\LCScale_11x11Grid\expectedoutput"
#sysag = "C:\iTree\HydroPlus\TestingFilesAndScript\TestCases\SpatialCoolAir\LCScale_12x12Grid\output"
#fpath = str(sysag  + '\\')

#change directory to file path with model inputs and where outputs are written
os.chdir(fpath)
path_out = fpath

#Obtain input from user via Command Prompt input
print('\n')
#cell_length is square array size, either 11 or 12 elements along a side
cell_length = input('Enter 11 for 11x11 or 12 for 12x12 to identify map array size (11 or 12):')
#r is Hourly or Daily air temperature metrics
#Note: Daily is typically used for HydroPlus test cases
r = input(
    'Enter H for Hourly or D for Daily Air Temperature Metric vs IC and TC Maps (H or D):')
    
#If r is H or h then set to H for error handling, and query about sampling duration
if r == 'H' or r == 'h':
    r = 'H'
    r2 = input(
        'Select the sampling duration of Air Temperature Metric Maps: Year (Y), Seasonal (S), Month (M), Day (D), or Hour (H):')
        
    #If r2 is H or h then set to H for error handling
    if r2 == 'H' or r2 == 'h':
        r2 = 'H'
    #If r2 is D or d then set to D for error handling
    if r2 == 'D' or r2 == 'd':
        r2 = 'D'
    #If r2 is M or m then set to M for error handling
    if r2 == 'M' or r2 == 'm':
        r2 = 'M'
    #If r2 is S or s then set to S for error handling
    if r2 == 'S' or r2 == 's':
        r2 = 'S'
    #If r2 is Y or y then set to Y for error handling
    if r2 == 'Y' or r2 == 'y':
        r2 = 'Y'        
#If r is D or d then set to D for error handling
if r == 'D' or r == 'd':
    r = 'D'
    r2 = None
   
r3 = None

#Use user given cell size to compute cell number
cell_length = int(cell_length)
cell_num = cell_length * cell_length
    
# if user wants to have hourly, daily, monthly yearly and seasonally plots, read the BlockGroup 1-121 files into matrix
if r == "H":
    # read multiple files
    a = "BlockGroup"
    c = ".csv"

    # to know the date and time
    #cte Refactor to error handle if BlockGroup1.csv does not exist
    timeseries = np.loadtxt(path_out + "BlockGroup1.csv", delimiter=',', skiprows=1, usecols=[0])
    timeseries = timeseries.astype(int)
    starttime = int(timeseries[0])
    endtime = int(timeseries[-1])
    print("start time", starttime, "end time", endtime)

    year = timeseries // 1000000
    month = (timeseries // 10000) % 10
    day = (timeseries // 100) % 100
    time = timeseries % 100

    uniq_year = np.unique(year)

    # create a matrix to save all the temperature data (each hour 121/144 data in a row)
    matrix = np.empty([len(timeseries), cell_num])
    for b in range(0, cell_num):
        file_name = a + str(b + 1) + c
        file_to_open = path_out + file_name
        tempT = np.loadtxt(file_to_open, delimiter=',', skiprows=1, usecols=[1])
        matrix[:, b] = tempT

    # create a matrix with time, with column Year, Month, Day, Time, 121/144 data)
    matrix_w_date = np.column_stack((year, month, day, time, matrix))

    print('Finished reading BlockGroup files')

# if user wants to have max or min plots, read the BlockGroupDaily 1-121/144 files into matrix2
if (r == "D"):
    # read multiple files
    a = "BlockGroupDaily"
    c = ".csv"

    # to know the date and time
    days = np.loadtxt(path_out + "BlockGroupDaily1.csv", delimiter=',', skiprows=1, usecols=[0])
    # if just 1 day, the syntax days[0] will cause error, so we separate it into 2 senario,=1 or >1
    if days.size == 1:
        days = np.array([days.astype(int)])
    if days.size > 1:
        days = np.array(days.astype(int))
    startdate = days[0]
    enddate = days[-1]
    print("Start date (YYYYMMDD):", startdate)
    print("End date (YYYYMMDD):", enddate)

    year = days // 10000
    month = (days // 100) % 10
    day = days % 100

    uniq_year = np.unique(year)

    # create a matrix2 to save all the temperature data (121/144 data in a row)
    # RawTairMax_K 121/144 values
    # RawTairAvg_K 121/144 values
    # ...
    matrix2 = np.empty([11, cell_num])
    # save 11 max or min dates
    maxmindate = np.empty([cell_num, 11])

    for b in range(0, cell_num):
        file_name = a + str(b + 1) + c
        file_to_open = path_out + file_name
        # read each BlockGroupDaily file and put it in TempT (with header)
        tempT = np.loadtxt(file_to_open, delimiter=',', skiprows=1, usecols=range(1, 12))

        for i in range(0, 11):
            if days.size == 1:
                col_name = tempT[i]
                tempatok = tempT[i]
                matrix2[i, b] = tempatok
                maxmindate[b, i] = days

            if days.size > 1:
                # get the Max of these columns and put them into rows in matrix2
                if (i == 0) | (i == 2) | (i == 3) | (i == 4) | (i == 5) | (i == 6) | (i == 7) | (i == 8):
                    tempatok = np.amax(tempT[:, i])
                    maxmindate[b, i] = days[np.argmax(tempT[:, i])]
                    # get the Min of these columns and put them into rows in matrix2
                elif (i == 1) | (i == 9) | (i == 10):
                    tempatok = np.amin(tempT[:, i])
                    maxmindate[b, i] = days[np.argmin(tempT[:, i])]
                matrix2[i, b] = tempatok
    # matrix2_w_date=np.column_stack((days,matrix2))
    maxmindate.astype('int')
    
    print('Finished reading BlockGroupDaily files')
    # txtpath=path_out+"/%s.csv" % ("matrix2")
    # np.savetxt(txtpath,matrix2,fmt="%f")
    
########################################
#             Functions                #
########################################
#create folders path to save the plots
def createfolder (foldername):
    outputpath = os.path.join(path_out, foldername)
    isdir = os.path.isdir(outputpath)
    if not isdir:
        os.makedirs(outputpath)
        print("Directory '% s' created" % foldername)
    else: print("Directory '% s' existed, saving files to this folder..." % foldername)  
    return outputpath
    
# to find the data for that hour in the matrix
def find_data(date,row, column):
    y = row * cell_length + column 
    y = int(y)
    if (r == 'H') & (r2 == 'H'):
        z = matrix[date,y]
    if (r == 'H') & (r2 == 'Y'):
        z = matrix_year[0, y]
    if (r == 'H') & (r2 == 'M'):    
        z = matrix_month[0, y]
    if (r == 'H') & (r2 == 'D'):
        z = matrix_day[0, y]
    if (r == 'H') & (r2 == 'S'):
        z = matrix_season[0, y]
    if (r == 'D'):
        z = matrix2[date,y]
    
    z.astype(float)
    return z
    
def f(x,y):
    return Z_m[x, y]
      
#plot xyz prepare and plot
def plot_3D(cell_length,xlabel,ylabel,zlabel,outputpath):
    if (cell_length == 11):
        Xarr = np.arange( 0, cell_length *10, 10)
        Yarr = np.arange( 0, cell_length *10, 10)
        X, Y = np.meshgrid(Xarr, Yarr)
        Z = f(Y//10,X//10)
        
    if (cell_length == 12):
        Xarr = (0, 1 ,10 ,20 ,30 ,40 ,50 ,60 ,70 ,80 ,90 ,100)
        Yarr = (0, 1 ,10 ,20 ,30 ,40 ,50 ,60 ,70 ,80 ,90 ,100)
        X, Y = np.meshgrid(Xarr, Yarr)
        Z = Z_m
        
    fig = plt.figure()
    ax = fig.add_subplot(projection='3d')
    surf = ax.plot_surface(X, Y, Z, cmap=cm.coolwarm, linewidth=0, antialiased=False)
    ax.invert_yaxis()
    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)
    ax.set_zlabel(zlabel)

    # Customize the z axis.
    ax.zaxis.set_major_locator(LinearLocator(10))
    ax.zaxis.set_major_formatter(FormatStrFormatter('%.02f'))
    fig.colorbar(surf, shrink=0.5, aspect=5)
    if (r == 'H'):
        if(r2 == 'H'):
            figpath = outputpath + "/%s.png" % (timeseries[i])
        if(r2 == 'Y'):
            figpath = outputpath + "/%s.png" % (uniq_year[j])
        if(r2 == 'M'): 
            figpath = outputpath + "/%s_%s.png" % (uniq_year[i], (j + 1))
        if (r2 == 'D'):
            figpath = outputpath + "/%s_%s_%s.png" % (uniq_year[i], j + 1, k + 1)
    if (r == 'H') & (r2 == 'S'):
        if (season == "summer"):
            figpath = outputpath + "/%s_%s.png" % (uniq_year[i], 'summer')
        if(season == "spring"):
            figpath = outputpath + "/%s_%s.png" % (uniq_year[i], 'spring')
        if(season == "autumn"):
            figpath = outputpath + "/%s_%s.png" % (uniq_year[i], 'autumn')
        if(season == "winter"):
            figpath = outputpath + "/%s_%s.png" % (uniq_year[i], 'winter')
    if (r == 'D'): 
        figpath = outputpath + "/%s.png" % (name[i])
    if (r3 == "Y"):
        ax.view_init(elev=30, azim=-60)
        plt.show()
    else:
        ax.view_init(elev=30, azim=-60)
        plt.savefig(figpath, format='png') 
        plt.close()
########################################
# Main function
# start creating plots for each hour
if (r == 'H') & (r2 == 'H'):
    hr_outputpath = createfolder(str("hourly"))

    # create plots for each hour in floder hourlyplot
    for i in range(0, len(timeseries)):
        # create a new matrix to save the 11*11 or 12*12 2d data
        Z_m = np.empty([cell_length, cell_length])
        for X in range(0, cell_length, 1):
            for Y in range(0, cell_length, 1):
                Z_m[X, Y] = find_data(i, int(X), int(Y))
   
        plot_3D(cell_length, 'TC','IC','Temperature', hr_outputpath)
        
    print('Saved plots for each hour')


# create plots for each year
if (r == 'H') & (r2 == 'Y'):
    y_outputpath = createfolder(str("yearly"))
    
    for j in range(0, len(uniq_year)):
        matrix_filter = matrix_w_date[matrix_w_date[:, 0] == uniq_year[j]]
        if len(matrix_filter) != 0:
            matrix_year = np.empty([1, cell_num])
            for k in range(0, cell_num):
                matrix_year[:, k] = np.mean(matrix_filter[:, k + 4])

        Z_m = np.empty([cell_length, cell_length])
        for X in range(0, cell_length, 1):
            for Y in range(0, cell_length, 1):
                Z_m[X,Y] = find_data(0, int(X), int(Y))
        plot_3D(cell_length, 'TC','IC','Temperature', y_outputpath)
    print('Saved plots for each year (averaged)')


# create plots for each month
if (r == 'H') & (r2 == 'M'):
    m_outputpath = createfolder(str("monthly"))

    for i in range(0, len(uniq_year)):
        for j in range(0, 12):
            matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (matrix_w_date[:, 1] == j + 1)]
            if len(matrix_filter) != 0:
                matrix_month = np.empty([1, cell_num])
                for k in range(0, cell_num):
                    matrix_month[:, k] = np.mean(matrix_filter[:, k + 4])

                # put the data we found in a new 11*11 matrix Z_m for plotting
                Z_m = np.empty([cell_length, cell_length])
                for X in range(0, cell_length, 1):
                    for Y in range(0, cell_length, 1):
                        Z_m[X,Y] = find_data(0, int(X), int(Y))
            
                plot_3D(cell_length, 'TC','IC','Temperature', m_outputpath)                   
    print('Saved plots for each month (averaged)')


# create plots for each day
if (r == 'H') & (r2 == 'D'):
    d_outputpath = createfolder(str("daily"))
    
    for i in range(0, len(uniq_year)):
        for j in range(0, 12):
            for k in range(0, 31):
                matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (matrix_w_date[:, 1] == j + 1) & (
                            matrix_w_date[:, 2] == k + 1)]
                if len(matrix_filter) != 0:
                    matrix_day = np.empty([1, cell_num])
                    for l in range(0, cell_num):
                        matrix_day[0, l] = np.mean(matrix_filter[:, l + 4])

                    # put the data we found in a new 11*11 matrix Z_m for plotting
                    Z_m = np.empty([cell_length, cell_length])
                    for X in range(0, cell_length, 1):
                        for Y in range(0, cell_length, 1):
                            Z_m[X,Y] = find_data(0, int(X), int(Y))
                    
                    plot_3D(cell_length, 'TC','IC','Temperature', d_outputpath)     

    print('Saved plots for each day (averaged)')

# create plots for each season
# create matrix to save averaged data for each year
if (r == 'H') & (r2 == 'S'):
    s_outputpath = createfolder(str("season"))
    
    # set up the filter to find the day that year and month
    for i in range(0, len(uniq_year)):
        matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (
                    (matrix_w_date[:, 1] == 3) | (matrix_w_date[:, 1] == 4) | (matrix_w_date[:, 1] == 5))]
        season = "spring"
        if len(matrix_filter) != 0:
            matrix_season = np.empty([1, cell_num])
            for k in range(0, cell_num):
                matrix_season[:, k] = np.mean(matrix_filter[:, k + 4])

            # put the data we found in a new 11*11 matrix Z_m for plotting
            Z_m = np.empty([cell_length, cell_length])
            for X in range(0, cell_length, 1):
                for Y in range(0, cell_length, 1):
                    Z_m[X,Y] = find_data(0, int(X), int(Y))
            plot_3D(cell_length, 'TC','IC','Temperature', s_outputpath) 
            print("Saved plots in output folder for '%s' spring" % uniq_year[i])

    for i in range(0, len(uniq_year)):
        matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (
                    (matrix_w_date[:, 1] == 6) | (matrix_w_date[:, 1] == 7) | (matrix_w_date[:, 1] == 8))]
        season = "summer"
        if len(matrix_filter) != 0:
            matrix_season = np.empty([1, cell_num])
            for k in range(0, cell_num):
                matrix_season[:, k] = np.mean(matrix_filter[:, k + 4])

            # put the data we found in a new 11*11 matrix Z_m for plotting
            Z_m = np.empty([cell_length, cell_length])
            for X in range(0, cell_length, 1):
                for Y in range(0, cell_length, 1):
                    Z_m[X,Y] = find_data(0, int(X), int(Y))
            plot_3D(cell_length, 'TC','IC','Temperature', s_outputpath) 
            print("Saved plots in output folder for '%s' summer" % uniq_year[i])

    for i in range(0, len(uniq_year)):
        matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (
                    (matrix_w_date[:, 1] == 9) | (matrix_w_date[:, 1] == 10) | (matrix_w_date[:, 1] == 11))]
        season = "autumn"
        if len(matrix_filter) != 0:
            matrix_season = np.empty([1, cell_num])
            for k in range(0, cell_num):
                matrix_season[:, k] = np.mean(matrix_filter[:, k + 4])

            # put the data we found in a new 11*11 matrix Z_m for plotting
            Z_m = np.empty([cell_length, cell_length])
            for X in range(0, cell_length, 1):
                for Y in range(0, cell_length, 1):
                    Z_m[X,Y] = find_data(0, int(X), int(Y))
            plot_3D(cell_length, 'TC','IC','Temperature', s_outputpath) 
            print("Saved plots in output folder for '%s' autumn" % uniq_year[i])

    for i in range(0, len(uniq_year)):
        matrix_filter = matrix_w_date[(matrix_w_date[:, 0] == uniq_year[i]) & (
                    (matrix_w_date[:, 1] == 12) | (matrix_w_date[:, 1] == 1) | (matrix_w_date[:, 1] == 2))]
        seanson = "winter"
        if len(matrix_filter) != 0:
            matrix_season = np.empty([1, cell_num])
            for k in range(0, cell_num):
                matrix_season[:, k] = np.mean(matrix_filter[:, k + 4])


            # a function for finding the data for each blockgroup in matrix
            Z_m = np.empty([cell_length, cell_length])
            for X in range(0, cell_length, 1):
                for Y in range(0, cell_length, 1):
                    Z_m[X,Y] = find_data(0, int(X), int(Y))
            plot_3D(cell_length, 'TC','IC','Temperature', s_outputpath) 
            print("Saved plots in output folder for '%s' winter" % uniq_year[i])

########################################
# plot for file BlockGroupDaily#.csv

if (r == "D"):
    atok_outputpath = createfolder(str("Daily_Metrics"))
    
    name = np.array(['Tair_Max_K', 'Tair_Min_K', 'Tair_Avg_K', 'Tdew_Avg_K', 'HeatIndex_NWS_Max_K',
                     'HeatIndex_NWS_Avg_K', 'HeatIndex_Lin_Max_K', 'HeatIndex_Lin_Avg_K',
                     'HeatIndex_NWS_AvgDaily_K', 'WindChill_NWS_Min_K', 'WindChill_NWS_Avg_K'])
    for i in range(0, 11):   
        Z_m = np.empty([cell_length, cell_length])
        for X in range(0, cell_length, 1):
            for Y in range(0, cell_length, 1):
                Z_m[X,Y] = find_data(i, int(X), int(Y))
        
        plot_3D(cell_length, 'TC','IC','Temperature', atok_outputpath) 
        
    print("Saved plots in output folder for Daily_Metrics")

    #########################################
    # generate csv files
    header = ['Block', 'Tair_Max_K', 'Tair_Min_K', 'Tair_Avg_K', 'Tdew_Avg_K', 'HeatIndex_NWS_Max_K',
                     'HeatIndex_NWS_Avg_K', 'HeatIndex_Lin_Max_K', 'HeatIndex_Lin_Avg_K',
                     'HeatIndex_Lin_AvgDaily_K', 'WindChill_NWS_Min_K', 'WindChill_NWS_Avg_K']
    header_v = np.array([np.arange(1, cell_num+1, dtype= int)])
    matrix_date = np.column_stack((header_v.T, maxmindate.astype('int')))
    matrix_date = np.vstack((header, matrix_date))

    info_headers = "a.Tair_Max_K – plot the maximum\r\nb.Tair_Min_K – plot the minimum\r\nc.Tair_Avg_K – plot the average\r\nd.Tdew_Avg_K – plot the average\r\ne.HeatIndex_NWS_Max_K – plot the maximum\r\nf.HeatIndex_NWS_Avg_K – plot the average\r\ng.HeatIndex_Lin_Max_K – plot the maximum\r\nh.HeatIndex_Lin_Avg_K – plot the average\r\ni.HeatIndex_Lin_AvgDaily_K – plot the average\r\nj.WindChill_NWS_Min_K – plot the minimum\r\nk.WindChill_NWS_Avg_K – plot the average\r\n"

    txtpath = atok_outputpath + "/%s.csv" % ("Daily_Metrics")
    np.savetxt(txtpath, matrix_date, fmt="%s", delimiter=',', header=info_headers, comments='')
    print("Daily_Metrics.csv file saved to Daily_Metrics folder")

    info_date = np.empty([11, 2])
    for i in range(0, 11):
        repeatedmax = 0
        datemax = 0
        uniquedate = np.unique(maxmindate[:, i])
        for j in range(0, len(uniquedate)):
            date_from_maxmin = maxmindate[:, i]
            # date_from_maxmin.astype('int')
            if np.count_nonzero(date_from_maxmin == uniquedate[j]) > repeatedmax:
                repeatedmax = np.count_nonzero(date_from_maxmin == uniquedate[j])
                datemax = uniquedate[j]
        info_date[i, 0] = datemax
        info_date[i, 1] = repeatedmax

    info = np.empty([11, 3])

    info[:, 0] = np.repeat(len(days), 11)
    info[:, 1] = info_date[:, 0]
    info[:, 2] = info_date[:, 1]
    MaxorMin = np.array(['Max', 'Min', 'Max', 'Max', 'Max', 'Max', 'Max', 'Max', 'Max', 'Min', 'Min'])
    header_v2 = np.array(
        ['Tair_Max_K.png', 'Tair_Min_K.png', 'Tair_Avg_K', 'Tdew_Avg_K.png', 'HeatIndex_NWS_Max_K.png', 'HeatIndex_NWS_Avg_K.png',
         'HeatIndex_Lin_Max_K.png', 'HeatIndex_Lin_Avg_K.png', 'HeatIndex_Lin_AvgDaily_K.png', 'WindChill_NWS_Min_K.png',
         'WindChill_NWS_Avg_K.png'])
    header2 = [["*Number of blocks reached the max/min T on that day", "", "", "", ""],
               ["filename", "UniqueDates", "YYYYMMDD", "Blocks", "Type"]]

    matrix_info = np.column_stack((header_v2.T, info.astype('int'), MaxorMin.T))
    matrix_info = np.vstack((header2, matrix_info))
    txtpath = atok_outputpath + "/%s.csv" % ('AirTempMetric-Info')
    np.savetxt(txtpath, matrix_info, fmt="%s", delimiter=',', comments='')
    print('AirTempMetric-Info.csv file saved to Daily_Metrics folder')

    flag = True
    while flag == True:
        r3 = input('Do you want to rotate a 3D plot and save it manually? (Y/N):')
        if r3 == "Y" or r3 == 'y':
            r3 = 'Y'
            
            flag = True
            r4 = input(
                '3D plot types available to rotate:\r\nMaximun Tair_Max_K (a)\r\nMinimun Tair_Min_K (b)\r\nMaverage Tair_Avg_K (c)\r\nAverage Tdew_Avg_K (d)\r\nMaximum HeatIndex_NWS_Max_K (e)\r\nAverage HeatIndex_NWS_Avg_K (f)\r\nMaximum HeatIndex_Lin_Max_K (g)\r\nAverage HeatIndex_Lin_Avg_K (h)\r\nAverage HeatIndex_Lin_AvgDaily_K (i)\r\nMinimum WindChill_NWS_Min_K (j)\r\nAverage WindChill_NWS_Avg_K (k)\r\nPlease choose a 3D plot type by letter (a to k):')
            if r4 == 'a' or  r4 == 'A':
                rowatok = 0
            if r4 == 'b' or  r4 == 'B':
                rowatok = 1
            if r4 == 'c' or  r4 == 'C':
                rowatok = 2
            if r4 == 'd' or  r4 == 'D':
                rowatok = 3
            if r4 == 'e' or  r4 == 'E':
                rowatok = 4
            if r4 == 'f' or  r4 == 'F':
                rowatok = 5
            if r4 == 'g' or  r4 == 'G':
                rowatok = 6
            if r4 == 'h' or  r4 == 'H':
                rowatok = 7
            if r4 == 'i' or  r4 == 'I':
                rowatok = 8
            if r4 == 'j' or  r4 == 'J':
                rowatok = 9
            if r4 == 'k' or  r4 == 'K':
                rowatok = 10

            Z_m = np.empty([cell_length, cell_length])
            for X in range(0, cell_length, 1):
                for Y in range(0, cell_length, 1):
                    Z_m[X,Y] = find_data(rowatok, int(X), int(Y))

            plot_3D(cell_length, 'TC','IC','Temperature', "NoPathNeeded") 

            plt.show()
            # print('Showed')
            # plt.close()
            # print("Saved plots for '% s'"% directory)
        else:
            print('\n')
            print("Simulation completed. Check output folder for results.")
            print("Thank you for using i-Tree tools to improve the world!")
            flag = False
