EX1_model2+3.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "mosel2_3.ipynb",
      "provenance": [],
      "authorship_tag": "ABX9TyPWAF+1YjsUFiCtw7JT3tBp",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/sonomarina/Mathmodelling_DTU/blob/main/EX1_model2%2B3.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "pypbVqJxnf0x"
      },
      "source": [
        "# DATA FILES DIRECTORY\n",
        "dirIn = 'path to Data folder for exercise 1'\n",
        "# IMPORTS\n",
        "import helpFunctions as hf\n",
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import imageio\n",
        "# DAY ONE MULTISPECTRAL IMAGE AND ANNOTATIONS\n",
        "multiIm, annotationIm = hf.loadMulti('multispectral_day01.mat','annotation_day01.png',\\\n",
        "dirIn)\n",
        "# ASSUME EQUAL PROBABILITY. CHANGE THIS TO ACCOMMODATE METHOD 3 PRIOR PROBABILITIES\n",
        "p = 0.5\n",
        "# EXTRACT PIXELS FOR FAT AND MEAT\n",
        "[fatPix, fatR, fatC] = hf.getPix(multiIm, annotationIm[:,:,1]);\n",
        "[meatPix, meatR, meatC] = hf.getPix(multiIm, annotationIm[:,:,2]);\n",
        "# MEANS FOR THE BANDS\n",
        "meatM = np.mean(meatPix,0)\n",
        "fatM = np.mean(fatPix, 0)\n",
        "covM = np.zeros((19,19))\n",
        "covF = np.zeros((19,19))\n",
        "# ARRAY OF PIXELS FOR EACH BAND\n",
        "m = 264196\n",
        "reshIm = np.zeros((m, 19))\n",
        "for i in range(0, 19):\n",
        "reshIm[:,i] = multiIm[:,:,i].flatten()\n",
        "# COMPUTE COVARIANCE\n",
        "for i in range(0, 19):\n",
        "for j in range(0, 19):\n",
        "covM[i, j] = 1/(m - 1)*np.sum((reshIm[:,i] - meatM[i])*(reshIm[:,j] - meatM[j]))\n",
        "covF[i, j] = 1/(m - 1)*np.sum((reshIm[:,i] - fatM[i])*(reshIm[:,j] - fatM[j]))\n",
        "cov = 1/(2*m - 2)*((m-1)*covM + (m-1)*covF)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Mp4zvQ5Won2t"
      },
      "source": [
        "# LDA FUNCTION\n",
        "\n",
        "def discFunc(x):\n",
        "  s_meat = np.dot(x,np.dot(np.linalg.inv(cov),meatM)) - \\\n",
        "1/2*np.dot(meatM, np.dot(np.linalg.inv(cov),meatM)) + np.log(p)\n",
        "  s_fat = np.dot(x,np.dot(np.linalg.inv(cov),fatM)) - \\\n",
        "1/2*np.dot(fatM, np.dot(np.linalg.inv(cov),fatM)) + np.log(p)\n",
        "  return s_meat, s_fat\n",
        "# CLASSIFY THE PIXELS\n",
        "s = np.zeros(264196)\n",
        "for i in range(0, 264196):\n",
        "  s_meat, s_fat = discFunc(reshIm[i,:])\n",
        "  if (s_fat >= s_meat):\n",
        "    s[i] = 1\n",
        "  else:\n",
        "      s[i] = 0\n",
        "s = np.reshape(s, (514, 514)).astype(np.int64)\n",
        "# IDENTIFIED CLASSIFICATIONS VS KNOWN CLASSIFICATIONS\n",
        "c = s[fatR, fatC] == annotationIm[fatR, fatC, 1]\n",
        "rate = len(c[c == True])/len(c)\n",
        "# IDENTIFY BACKGROUND\n",
        "bPix = np.where(np.logical_and(annotationIm[:,:,0] == 0, annotationIm[:,:,1] == 0, \\\n",
        "annotationIm[:,:,2] == 0))\n",
        "# REMOVE BACKGROUND\n",
        "s[bPix] = 0\n",
        "# CLASSIFIED IMAGE PLOT\n",
        "[fatPix, fatR, fatC] = hf.getPix(multiIm, s)\n",
        "pixId = np.stack((fatR, fatC), axis=1)\n",
        "imRGB = imageio.imread(dirIn + 'color_day01.png')\n",
        "rgbOut = hf.setImagePix(imRGB, pixId)\n",
        "plt.imshow(rgbOut)\n",
        "plt.title('Classified image of the salami on day 1')\n",
        "plt.xlim(0,514)\n",
        "plt.ylim(514,0)\n",
        "plt.savefig('salami_day1_model_2', dpi = 300)\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}