# ==============================================================================
# Script name      : warehouse_manipulation.r
# Student name     : GOBINATH GOTHANDAPANI
# Student number   : 2538599
# ==============================================================================

# ##############################################################################
# Part A: Generate Warehouse Data
# ##############################################################################
# ---- A1. Set the random seed ------------------------------------------------
# Makes the "random" numbers the same every time the script runs
set.seed(123)

# ---- A2. Define the simulation parameters -----------------------------------
# Settings stored in one place so they are easy to change
number_of_days       <- 28     # 28 days
number_of_zones      <- 5      # 5 zones
number_of_categories <- 6      # 6 product categories
base_units           <- 100    # starting units
time_per_unit        <- 0.85   # minutes per unit

# Total number of rows = every combination of day, zone and category
total_rows <- number_of_days * number_of_zones * number_of_categories   # 28 x 5 x 6 = 840

# ---- A3. Create an empty dataframe with the correct columns -----------------
# An empty table of the right size (filled with 0s); the loops fill it in later
warehouse_data <- data.frame(
  Day            = integer(total_rows),
  Zone           = integer(total_rows),
  Category       = integer(total_rows),
  UnitsProcessed = integer(total_rows),
  ProcessingTime = numeric(total_rows)
)

# ---- A4. Populate the dataframe using nested loops --------------------------
# Three loops (day > zone > category) create one row for every combination
row_counter <- 1   # Tracks which row of the dataframe is being filled

for (day_id in 1:number_of_days) {                    # Loop 1: each day
  for (zone_id in 1:number_of_zones) {                # Loop 2: each zone
    for (category_id in 1:number_of_categories) {     # Loop 3: each category
      
      # RandomVariation: a random integer between -12 and +12 (inclusive)
      random_variation <- sample(-12:12, size = 1)
      
      # UnitsProcessed = BaseUnits + (Zone * 8) + (Category * 5) + RandomVariation
      units_processed <- base_units + (zone_id * 8) + (category_id * 5) + random_variation
      
      # ProcessingTime = UnitsProcessed * TimePerUnit
      processing_time <- units_processed * time_per_unit
      
      # Store the values in the current row of the dataframe
      warehouse_data$Day[row_counter]            <- day_id
      warehouse_data$Zone[row_counter]           <- zone_id
      warehouse_data$Category[row_counter]       <- category_id
      warehouse_data$UnitsProcessed[row_counter] <- units_processed
      warehouse_data$ProcessingTime[row_counter] <- processing_time
      
      # Move on to the next row
      row_counter <- row_counter + 1
    }
  }
}

# ---- A5. Display the first 10 rows ------------------------------------------
# head(data, 10) = show the first 10 rows
cat("First 10 rows of warehouse_data:\n")
print(head(warehouse_data, 10))

# ---- A6. Display the total number of rows and columns -----------------------
# nrow() = number of rows, ncol() = number of columns
cat("\nTotal number of rows   :", nrow(warehouse_data), "\n")
cat("Total number of columns:", ncol(warehouse_data), "\n")


# ##############################################################################
# Part B: Total Units Processed per Zone
# ##############################################################################
# ---- B1. Create the zone_totals dataframe -----------------------------------
# Results table: one row per zone
zone_totals <- data.frame(
  Zone       = 1:number_of_zones,
  TotalUnits = numeric(number_of_zones)
)

# ---- B2. Calculate each zone's total using loops ----------------------------
# For each zone, check every record and add up the units that belong to it
for (zone_id in 1:number_of_zones) {          # Outer loop: each zone
  
  running_total <- 0                          # Reset the total for this zone
  
  for (row_id in 1:total_rows) {              # Inner loop: every record
    
    # Only add the units when the record belongs to the current zone
    if (warehouse_data$Zone[row_id] == zone_id) {
      running_total <- running_total + warehouse_data$UnitsProcessed[row_id]
    }
  }
  
  # Save the finished total for this zone
  zone_totals$TotalUnits[zone_id] <- running_total
}

# ---- B3. Display the completed dataframe ------------------------------------
cat("Total units processed per zone:\n")
print(zone_totals)


# ##############################################################################
# Part C: Average Processing Time per Product Category
# ##############################################################################
# ---- C1. Create the category_average_time dataframe -------------------------
# Results table: one row per category
category_average_time <- data.frame(
  Category              = 1:number_of_categories,
  AverageProcessingTime = numeric(number_of_categories)
)

# ---- C2. Calculate each category's average using loops ----------------------
# Average = total time / number of records, so we add up both inside the loop
for (category_id in 1:number_of_categories) {   # Outer loop: each category
  
  total_processing_time <- 0   # Running total of processing time (minutes)
  record_count          <- 0   # Number of records found for this category
  
  for (row_id in 1:total_rows) {                # Inner loop: every record
    
    # Only use records that belong to the current category
    if (warehouse_data$Category[row_id] == category_id) {
      total_processing_time <- total_processing_time + warehouse_data$ProcessingTime[row_id]
      record_count          <- record_count + 1
    }
  }
  
  # Average = total processing time / number of records (rounded to 2 d.p.)
  category_average_time$AverageProcessingTime[category_id] <-
    round(total_processing_time / record_count, 2)
}

# ---- C3. Display the completed dataframe ------------------------------------
cat("Average processing time (minutes) per product category:\n")
print(category_average_time)


# ##############################################################################
# Part D: Identify Stock Processing Anomalies
# ##############################################################################
# Anomaly = units more than 15% above OR below the average of its day and zone
# Example: average 100 -> normal range is 85 to 115

# ---- D1. Calculate the average units for every day/zone combination ---------
# A matrix is used as a lookup table:
#   rows = days (1 to 28), columns = zones (1 to 5)
# daily_zone_average[3, 2] = average units for day 3, zone 2.
daily_zone_average <- matrix(0, nrow = number_of_days, ncol = number_of_zones)

for (day_id in 1:number_of_days) {               # Loop 1: each day
  for (zone_id in 1:number_of_zones) {           # Loop 2: each zone
    
    running_total <- 0   # Total units for this day and zone
    record_count  <- 0   # Number of records (should be 6, one per category)
    
    for (row_id in 1:total_rows) {               # Loop 3: every record
      
      # Only use records matching BOTH the current day and the current zone
      if (warehouse_data$Day[row_id] == day_id && warehouse_data$Zone[row_id] == zone_id) {
        running_total <- running_total + warehouse_data$UnitsProcessed[row_id]
        record_count  <- record_count + 1
      }
    }
    
    # Store the average for this day/zone combination
    daily_zone_average[day_id, zone_id] <- running_total / record_count
  }
}

# ---- D2. Create the processing_anomalies dataframe --------------------------
# Results table: each record, its group average and a TRUE/FALSE anomaly flag
processing_anomalies <- data.frame(
  Day              = integer(total_rows),
  Zone             = integer(total_rows),
  Category         = integer(total_rows),
  UnitsProcessed   = integer(total_rows),
  DailyZoneAverage = numeric(total_rows),
  Anomaly          = logical(total_rows)
)

# ---- D3. Fill the dataframe and flag the anomalies using a loop -------------
for (row_id in 1:total_rows) {
  
  # Read the details of the current record from warehouse_data
  current_day      <- warehouse_data$Day[row_id]
  current_zone     <- warehouse_data$Zone[row_id]
  current_category <- warehouse_data$Category[row_id]
  current_units    <- warehouse_data$UnitsProcessed[row_id]
  
  # Look up the average for this record's day and zone
  current_average <- daily_zone_average[current_day, current_zone]
  
  # Anomaly rule: more than 15% above OR more than 15% below the average
  # ( || means OR )
  if (current_units > current_average * 1.15 || current_units < current_average * 0.85) {
    is_anomaly <- TRUE
  } else {
    is_anomaly <- FALSE
  }
  
  # Store everything in the processing_anomalies dataframe
  processing_anomalies$Day[row_id]              <- current_day
  processing_anomalies$Zone[row_id]             <- current_zone
  processing_anomalies$Category[row_id]         <- current_category
  processing_anomalies$UnitsProcessed[row_id]   <- current_units
  processing_anomalies$DailyZoneAverage[row_id] <- round(current_average, 2)  # Rounded for display only
  processing_anomalies$Anomaly[row_id]          <- is_anomaly
}

# ---- D4. Count how many anomalies were found (using a loop) -----------------
# Add 1 to the counter each time Anomaly is TRUE
anomaly_count <- 0
for (row_id in 1:total_rows) {
  if (processing_anomalies$Anomaly[row_id] == TRUE) {
    anomaly_count <- anomaly_count + 1
  }
}
cat("Number of anomalies detected:", anomaly_count, "out of", total_rows, "records\n\n")

# ---- D5. Display the first 10 rows where Anomaly is TRUE --------------------
# Keep only the rows where Anomaly is TRUE, then show the first 10
anomalies_only <- processing_anomalies[processing_anomalies$Anomaly == TRUE, ]

cat("First 10 anomalous records:\n")
print(head(anomalies_only, 10), row.names = FALSE)


# ##############################################################################
# Part E: Identify the Most Productive Day for Each Zone
# ##############################################################################
# ---- E1. Create the peak_zone_days dataframe --------------------------------
# Results table: one row per zone
peak_zone_days <- data.frame(
  Zone       = 1:number_of_zones,
  Day        = integer(number_of_zones),
  TotalUnits = numeric(number_of_zones)
)

# ---- E2. Find each zone's most productive day using loops -------------------
# For each zone, add up each day's units and remember the highest total
for (zone_id in 1:number_of_zones) {            # Loop 1: each zone
  
  best_day   <- 0   # Day with the highest total found so far
  best_total <- 0   # Highest daily total found so far
  
  for (day_id in 1:number_of_days) {            # Loop 2: each day
    
    day_total <- 0                              # Total units for this zone on this day
    
    for (row_id in 1:total_rows) {              # Loop 3: every record
      
      # Add the units only when the record matches the current zone AND day
      if (warehouse_data$Zone[row_id] == zone_id && warehouse_data$Day[row_id] == day_id) {
        day_total <- day_total + warehouse_data$UnitsProcessed[row_id]
      }
    }
    
    # If this day beats the best so far, remember it.
    # Strict ">" means that if two days tie, the earlier day is kept.
    if (day_total > best_total) {
      best_total <- day_total
      best_day   <- day_id
    }
  }
  
  # Save the peak day for this zone
  peak_zone_days$Day[zone_id]        <- best_day
  peak_zone_days$TotalUnits[zone_id] <- best_total
}

# ---- E3. Display the completed dataframe ------------------------------------
cat("Most productive day for each zone:\n")
print(peak_zone_days)


# ##############################################################################
# Optional Sanity Check (Loop-Based)
# ##############################################################################
# Add up all units in two ways. If both totals match, the zone totals are correct.

# Grand total taken directly from the raw data
grand_total_raw <- 0
for (row_id in 1:total_rows) {
  grand_total_raw <- grand_total_raw + warehouse_data$UnitsProcessed[row_id]
}

# Grand total built from the zone_totals dataframe
grand_total_zones <- 0
for (zone_id in 1:number_of_zones) {
  grand_total_zones <- grand_total_zones + zone_totals$TotalUnits[zone_id]
}

cat("Grand total from warehouse_data:", grand_total_raw, "\n")
cat("Grand total from zone_totals   :", grand_total_zones, "\n")

if (grand_total_raw == grand_total_zones) {
  cat("Check passed: the zone totals match the raw data.\n")
} else {
  cat("Check FAILED: the totals do not match - review Part B.\n")
}


# ##############################################################################
# Validation Checks
# ##############################################################################
# ==============================================================================
# VALIDATION CHECKS
# Uses loops and if statements only: no apply family, no user-defined functions.
# ==============================================================================
# Each test prints PASS or FAIL. We count how many tests run and how many pass.
tests_run    <- 0
tests_passed <- 0

# ---- T1. Dimensions ----------------------------------------------------------
# Check: the table has 840 rows and 5 columns
tests_run <- tests_run + 1
if (nrow(warehouse_data) == 840 && ncol(warehouse_data) == 5) {
  tests_passed <- tests_passed + 1
  cat("T1 PASS: warehouse_data has 840 rows and 5 columns\n")
} else {
  cat("T1 FAIL: unexpected dimensions\n")
}

# ---- T2. Every Day-Zone-Category combination appears once, in order ---------
# Check: every combination is in its expected row position.
# Rows were written in nested-loop order, so the row for (day, zone, category)
# must be at position (day-1)*30 + (zone-1)*6 + category.
mismatches <- 0
for (day_id in 1:number_of_days) {
  for (zone_id in 1:number_of_zones) {
    for (category_id in 1:number_of_categories) {
      expected_row <- (day_id - 1) * number_of_zones * number_of_categories +
        (zone_id - 1) * number_of_categories + category_id
      if (warehouse_data$Day[expected_row]      != day_id ||
          warehouse_data$Zone[expected_row]     != zone_id ||
          warehouse_data$Category[expected_row] != category_id) {
        mismatches <- mismatches + 1
      }
    }
  }
}
tests_run <- tests_run + 1
if (mismatches == 0) {
  tests_passed <- tests_passed + 1
  cat("T2 PASS: all 840 combinations are present exactly once\n")
} else {
  cat("T2 FAIL:", mismatches, "rows are out of place\n")
}

# ---- T3. Each day-zone group contains exactly 6 records ---------------------
# Check: every day/zone pair has 6 rows (one per category)
bad_groups <- 0
for (day_id in 1:number_of_days) {
  for (zone_id in 1:number_of_zones) {
    group_count <- 0
    for (row_id in 1:total_rows) {
      if (warehouse_data$Day[row_id] == day_id && warehouse_data$Zone[row_id] == zone_id) {
        group_count <- group_count + 1
      }
    }
    if (group_count != number_of_categories) {
      bad_groups <- bad_groups + 1
    }
  }
}
tests_run <- tests_run + 1
if (bad_groups == 0) {
  tests_passed <- tests_passed + 1
  cat("T3 PASS: all 140 day-zone groups contain 6 records\n")
} else {
  cat("T3 FAIL:", bad_groups, "groups have the wrong size\n")
}

# ---- T4. Formula, range and missing-value check for every row ---------------
# Check: no missing values, variation is between -12 and +12,
# ProcessingTime = Units * 0.85, and units stay between 101 and 182
bad_rows     <- 0
lowest_units <- warehouse_data$UnitsProcessed[1]
highest_units <- warehouse_data$UnitsProcessed[1]
for (row_id in 1:total_rows) {
  units_here <- warehouse_data$UnitsProcessed[row_id]
  time_here  <- warehouse_data$ProcessingTime[row_id]
  if (is.na(units_here) || is.na(time_here)) {      # is.na() = is the value missing?
    bad_rows <- bad_rows + 1
  } else {
    fixed_part <- base_units + warehouse_data$Zone[row_id] * 8 + warehouse_data$Category[row_id] * 5
    variation  <- units_here - fixed_part           # what is left must be the random part
    if (variation < -12 || variation > 12 || abs(time_here - units_here * time_per_unit) > 0.000000001) {
      bad_rows <- bad_rows + 1
    }
    if (units_here < lowest_units)  { lowest_units  <- units_here }
    if (units_here > highest_units) { highest_units <- units_here }
  }
}
tests_run <- tests_run + 1
if (bad_rows == 0 && lowest_units >= 101 && highest_units <= 182) {
  tests_passed <- tests_passed + 1
  cat("T4 PASS: formula holds in every row; units range from", lowest_units, "to", highest_units, "(allowed 101 to 182)\n")
} else {
  cat("T4 FAIL:", bad_rows, "bad rows; range", lowest_units, "to", highest_units, "\n")
}

# ---- T5. Zone totals recalculated by a different method (single pass) ------
# Check: recalculate zone totals with one pass through the data and compare
zone_check <- numeric(number_of_zones)
for (row_id in 1:total_rows) {
  zone_here <- warehouse_data$Zone[row_id]
  zone_check[zone_here] <- zone_check[zone_here] + warehouse_data$UnitsProcessed[row_id]
}
zone_errors <- 0
for (zone_id in 1:number_of_zones) {
  if (zone_check[zone_id] != zone_totals$TotalUnits[zone_id]) { zone_errors <- zone_errors + 1 }
}
tests_run <- tests_run + 1
if (zone_errors == 0) {
  tests_passed <- tests_passed + 1
  cat("T5 PASS: zone_totals match an independent single-pass recalculation\n")
} else {
  cat("T5 FAIL:", zone_errors, "zone totals differ\n")
}

# ---- T6. Category averages recalculated by a different method ---------------
# Check: recalculate category averages in one pass; each category has 140 records
category_time_check  <- numeric(number_of_categories)
category_count_check <- numeric(number_of_categories)
for (row_id in 1:total_rows) {
  category_here <- warehouse_data$Category[row_id]
  category_time_check[category_here]  <- category_time_check[category_here] + warehouse_data$ProcessingTime[row_id]
  category_count_check[category_here] <- category_count_check[category_here] + 1
}
category_errors <- 0
for (category_id in 1:number_of_categories) {
  recalculated <- round(category_time_check[category_id] / category_count_check[category_id], 2)
  if (recalculated != category_average_time$AverageProcessingTime[category_id] ||
      category_count_check[category_id] != 140) {
    category_errors <- category_errors + 1
  }
}
tests_run <- tests_run + 1
if (category_errors == 0) {
  tests_passed <- tests_passed + 1
  cat("T6 PASS: category averages match; each category has 140 records\n")
} else {
  cat("T6 FAIL:", category_errors, "categories differ\n")
}

# ---- T7. Results compared with analytic expected values ---------------------
# Check: results are close to what the maths predicts (within 3 SD = very likely range)
# Without noise a zone total is 28 * (705 + 48 * zone). Noise has SD about 7.2
# per record, so a zone total (168 records) should be within about 280 (3 SD).
# A category mean time is 0.85 * (124 + 5 * category), within about 1.6 (3 SD).
expected_errors <- 0
cat("\nT7 expected versus observed zone totals:\n")
for (zone_id in 1:number_of_zones) {
  expected_total <- 28 * (705 + 48 * zone_id)
  cat("  Zone", zone_id, "expected", expected_total, "observed", zone_totals$TotalUnits[zone_id], "\n")
  if (abs(zone_totals$TotalUnits[zone_id] - expected_total) > 280) { expected_errors <- expected_errors + 1 }
}
cat("T7 expected versus observed category average times:\n")
for (category_id in 1:number_of_categories) {
  expected_time <- 0.85 * (124 + 5 * category_id)
  cat("  Category", category_id, "expected", round(expected_time, 2), "observed", category_average_time$AverageProcessingTime[category_id], "\n")
  if (abs(category_average_time$AverageProcessingTime[category_id] - expected_time) > 1.6) { expected_errors <- expected_errors + 1 }
}
tests_run <- tests_run + 1
if (expected_errors == 0) {
  tests_passed <- tests_passed + 1
  cat("T7 PASS: all results lie within 3 SD of their expected values\n")
} else {
  cat("T7 FAIL:", expected_errors, "results outside tolerance (investigate; about 0.3% chance each by luck)\n")
}

# ---- T8. Peak days recalculated from a day-by-zone totals matrix ------------
# Check: build a day-by-zone totals matrix, then find each zone's best day again
day_zone_total <- matrix(0, nrow = number_of_days, ncol = number_of_zones)
for (row_id in 1:total_rows) {
  day_here  <- warehouse_data$Day[row_id]
  zone_here <- warehouse_data$Zone[row_id]
  day_zone_total[day_here, zone_here] <- day_zone_total[day_here, zone_here] + warehouse_data$UnitsProcessed[row_id]
}
peak_errors <- 0
for (zone_id in 1:number_of_zones) {
  check_day   <- 1
  check_total <- day_zone_total[1, zone_id]
  for (day_id in 2:number_of_days) {
    if (day_zone_total[day_id, zone_id] > check_total) {
      check_total <- day_zone_total[day_id, zone_id]
      check_day   <- day_id
    }
  }
  if (check_day != peak_zone_days$Day[zone_id] || check_total != peak_zone_days$TotalUnits[zone_id]) {
    peak_errors <- peak_errors + 1
  }
}
tests_run <- tests_run + 1
if (peak_errors == 0) {
  tests_passed <- tests_passed + 1
  cat("T8 PASS: peak_zone_days match the independent recalculation\n")
} else {
  cat("T8 FAIL:", peak_errors, "zones differ\n")
}

# ---- T9. Anomaly flags recalculated from the same day-by-zone matrix --------
# Check: recalculate every anomaly flag using the matrix from T8
flag_errors <- 0
for (row_id in 1:total_rows) {
  group_average <- day_zone_total[processing_anomalies$Day[row_id], processing_anomalies$Zone[row_id]] / 6
  units_here    <- processing_anomalies$UnitsProcessed[row_id]
  if (units_here > group_average * 1.15 || units_here < group_average * 0.85) {
    expected_flag <- TRUE
  } else {
    expected_flag <- FALSE
  }
  if (expected_flag != processing_anomalies$Anomaly[row_id]) { flag_errors <- flag_errors + 1 }
}
tests_run <- tests_run + 1
if (flag_errors == 0) {
  tests_passed <- tests_passed + 1
  cat("T9 PASS: every anomaly flag matches the recalculation\n")
} else {
  cat("T9 FAIL:", flag_errors, "flags differ\n")
}

# ---- T10. Small test case with a known outcome ------------------------------
# Check: a tiny example where we know the answer.
# Units 120, 80 and four values of 100 have average 100 (bounds 115 and 85),
# so exactly two records (120 and 80) must be flagged.
test_units <- c(120, 80, 100, 100, 100, 100)
test_total <- 0
for (item_id in 1:6) { test_total <- test_total + test_units[item_id] }
test_average <- test_total / 6
test_flagged <- 0
for (item_id in 1:6) {
  if (test_units[item_id] > test_average * 1.15 || test_units[item_id] < test_average * 0.85) {
    test_flagged <- test_flagged + 1
  }
}
tests_run <- tests_run + 1
if (test_flagged == 2) {
  tests_passed <- tests_passed + 1
  cat("T10 PASS: known test case flags exactly 2 records\n")
} else {
  cat("T10 FAIL: flagged", test_flagged, "records, expected 2\n")
}

# ---- Known limitation (boundary behaviour) ----------------------------------
# Computers store decimals approximately, so 100 * 1.15 might not equal 115 exactly.
# This only matters for a value sitting exactly on the 15% limit.
cat("\nBoundary note: is 100 * 1.15 exactly equal to 115?", 100 * 1.15 == 115, "\n")

# Final result: ideally "10 of 10"
cat("\nChecks passed:", tests_passed, "of", tests_run, "\n")