
/**
 * EmissionView.java
 * Used to process emission data. Calculates quantum yield, radiative rate constants,
 * optical density corrected emission spectra
 *
 * Created: Dec 30, 2004
 *
 * @author Nathan Stevens
 * @version 0.4 6/25/07
 */

import java.util.*;
import java.io.File;

public class EmissionView {
  private EmissionViewDialog emissionViewDialog = null;
  
  // main constructor
  public EmissionView(EmissionViewDialog emissionViewDialog) {
    this.emissionViewDialog = emissionViewDialog;
  }
  
  // Method to do quantum yield calculation
  public float[][] quantumYieldCalculation(Vector emData, Vector abData, 
  float[][] parameters) {
    int emSize = emData.size();
    int abSize = abData.size();
    
    float[][] results = new float[emSize][2];
    float ref_od = 1;
    float ref_area = 0;
    float ref_ri = (float)Math.pow(parameters[0][1], 2);  // squared refractive index of reference
    
    for(int i = 0; i < emSize; i++) {
      SpectraData emSD = (SpectraData)emData.elementAt(i);
      
      float sa_od = 1;
      if((i + 1) <= abSize) {
        SpectraData abSD = (SpectraData)abData.elementAt(i);
        sa_od = getOpticalDensity(abSD, parameters[i][2], parameters[i][5]);  // get the optcial density
      }
      
      // calculate quantum yield calculation
      if(i != 0) {
        float area  = getAreaUnderCurve(emSD, parameters[i][3], parameters[i][4]);
        float sa_ri = (float)Math.pow(parameters[i][1], 2);  // squared refractive index of sample
        
        float qy = parameters[0][0]*((area/ref_area)*(sa_ri/ref_ri)*(ref_od/sa_od));
        results[i][0] = qy;
        results[i][1] = area;
        
        /* debug code
        System.out.println("Standard Quantum Yield " + parameters[0][0]);
        System.out.println("Sample RI " + sa_ri + " Reference RI " + ref_ri);
        System.out.println("Sample OD " + sa_od + " Reference OD " + ref_od);
        System.out.println("Sample " +  emSD.getName() + " area under the curve is " + area);
        System.out.println("Sample " +  emSD.getName() + " Quantum Yield " + qy); //*/
      }
      else { // set the optical density for the reference
        ref_od = sa_od;
        ref_area = getAreaUnderCurve(emSD, parameters[0][3], parameters[0][4]);
        results[0][0] = 0.0f;
        results[0][1] = ref_area;
      }
    }
    
    return results;
  }
   
  // Method to do rate constant calculations
  public float[][] rateConstantCalculation(float[][] parameters) {
    float[][] results = new float[parameters.length][2];
    
    for(int i = 0; i < parameters.length; i++) {
      float qy = parameters[i][0];
      float lifetime = parameters[i][1];
      if(qy != 0.0f) {
        float kr = qy/lifetime;
        float knr = (1-qy)/lifetime;
        
        results[i][0] = kr;
        results[i][1] = knr;
        
        /*debug code
        System.out.println("Q = "  +  qy);
        System.out.println("Lifetime = "  + lifetime);
        System.out.println("kr = "  +  kr);
        System.out.println("knr = "  +  knr);*/
      }
      else {
        results[i][0] = 0;
        results[i][1] = 0;
      }
    }
    
    return results;
  }
  
  // Method to do optical density normalized emission
  public SpectraData[] odNormalizationCalculation(Vector emData, Vector abData, 
  float[][] parameters) {
    int size = emData.size();
    SpectraData[] results = new SpectraData[size];
    
    // used to normalize the data to one
    float max_int = 0.0f;
    float max_rint = 0.0f;
    
    for(int i = 0; i < size; i++) {
      SpectraData emSD = (SpectraData)emData.elementAt(i);
      SpectraData abSD = (SpectraData)abData.elementAt(i);
      
      String[] header = (String[])emSD.getHeader().clone();
      String[] dmeta = (String[])emSD.getMetaInfo().clone();
      String type = emSD.getType();
      
      float[][] data = emSD.getData();
      float[][] rdata = emSD.getReferenceData();
      float[][] ndata = new float[data.length][2];
      float[][] rndata = null;
      
      // set the reference data array if it is not null
      if(rdata != null) {
        rndata = new float[data.length][2];
      }
      
      // get the normalization factor
      float exWL = parameters[i][0];
      float nWL = parameters[i][1];
      
      float norm_factor = getOpticalDensity(abSD, exWL, nWL);
      
      for(int j = 0; j < data.length; j++) {
        ndata[j][0] = data[j][0];
        ndata[j][1] = data[j][1]/norm_factor;
        
        if(ndata[j][1] > max_int) {
          max_int = ndata[j][1];
        }
        
        if(rdata != null) {
           rndata[j][0] = rdata[j][0];
           rndata[j][1] = rdata[j][1]/norm_factor;
           
           if(rndata[j][1] > max_rint) {
             max_rint = rndata[j][1];
           }
        }
      }
      
      String name = emSD.getName() + "-ODN";
      results[i] = new SpectraData(name, type, header, dmeta, ndata, rndata);
    }
    
    // now normalize data to one
    for(int i = 0; i < results.length; i++) {
      SpectraData sd = results[i];
      
      float[][] data = sd.getData();
      float[][] rdata = sd.getReferenceData();
      float[][] ndata = new float[data.length][2];
      float[][] rndata = null;
      
      if(rdata != null) {
        rndata = new float[data.length][2];
      }
      
      for(int j = 0; j < data.length; j++) {
        ndata[j][0] = data[j][0];
        ndata[j][1] = data[j][1]/max_int;
        
        if(rdata != null) {
           rndata[j][0] = rdata[j][0];
           rndata[j][1] = rdata[j][1]/max_rint;
        }
      }
      
      sd.setData(ndata);
      sd.setReferenceData(rndata);
    }
    
    //System.out.println("Doing OD normalization calculation");  // debug code
    return results;
  }
  
  // Method to calculate the integrate intensity real and relative
  public float[][] getRelativeIntensity(SpectraData[] sda, float minWL, float maxWL) {
    float[][] relative_intensity = new float[sda.length][2];
    float area = 0.0f;
    float max_area = 0.0f;
    
    for(int i = 0; i < sda.length; i++) {
      area = getAreaUnderCurve(sda[i], minWL, maxWL);
      relative_intensity[i][0] = area;
      
      if(area > max_area) {
        max_area = area;
      }
    }
    
    // normalize the intergrated intensities to one now
    for(int i = 0; i < relative_intensity.length; i++) {
      relative_intensity[i][1] = relative_intensity[i][0]/max_area;
      
      /* debug code
      System.out.println(sda[i].getName() + " " + relative_intensity[i][0] + " : " + relative_intensity[i][1]); //*/
    }
    
    return relative_intensity;
  }
  
  // Method to return the relative intensity at certain wavelength
  public float[][] getRelativeWavelengthIntensity(SpectraData[] sda, float[] wls) {
    float[][] relative_intensity = new float[sda.length][wls.length];
    float intensity = 0.0f;
    float[] max_intensity = new float[wls.length];
    
    for(int i = 0; i < sda.length; i++) {
      for(int j = 0; j < wls.length; j++) {
        intensity = getOpticalDensity(sda[i], wls[j], 0.0f);
        relative_intensity[i][j] = intensity;
        
        if(intensity > max_intensity[j]) { // set the maximum intensity at this wavelength
          max_intensity[j] = intensity;
        }
      }
    }
    
    // nornamizled the results to 1 now
    for(int i = 0; i < sda.length; i++) {
      for(int j = 0; j < wls.length; j++) {
        relative_intensity[i][j] = relative_intensity[i][j]/max_intensity[j];
      }
    }
    
    return relative_intensity;
  }
  
  // function to divide by the minimum
  public float[][] divideByMinimum(float[][] data, String type) {
    if(type.equals("area")) {
      float min_area = data[0][1]; // set it to the firt data in the array
      float area = 0.0f;
      
      // find minimum
      for(int i = 0; i < data.length; i++) {
        area = data[i][1];
        
        if(area < min_area) {
          min_area = area;
        }
      }
      
      for(int i = 0; i < data.length; i++) {
        data[i][1]  /= min_area;
      }
    }
    else if(type.equals("wavelength")) {
      float intensity = 0.0f;
      float[] min_intensity = new float[data[0].length];
      
      // fill the minium intensity array
      for(int i = 0; i < data[0].length; i++) {
        min_intensity[i] = data[0][i];
      }
      
      // find minimum
      for(int i = 0; i < data.length; i++) {
        for(int j = 0; j < data[0].length; j++) {
          intensity = data[i][j];
          
          if(intensity < min_intensity[j]) { // set the maximum intensity at this wavelength
            min_intensity[j] = intensity;
          }
        }
      }
      
      // nornamizled the results to 1 now
      for(int i = 0; i < data.length; i++) {
        for(int j = 0; j < data[0].length; j++) {
          data[i][j] /=min_intensity[j];
        }
      }
    }
    else {
      System.err.println("Unkown data type.  Unable to divide by minimum");
    }
    
    return data;
  }
  
  // Method to return a normalized sensor
  public float[][] normalizeSensorData(float[][] data, String type) {
    if(type.equals("area")) {
      for(int i = 0; i < data.length; i++) {
        float old_value = data[i][1];
        data[i][1] /= old_value;  // makes this 1
        if(i + 1 != data.length) {
          data[i+1][1] = data[i+1][1]/old_value;  // normailze this
        }
        i++; // increment one more time
      }
    }
    else if(type.equals("wavelength")) { // must code
      // nornamizled the results to 1 now
      for(int i = 0; i < data.length; i++) {
        for(int j = 0; j < data[0].length; j++) {
          float old_value = data[i][j];
          data[i][j] /= old_value;
          if(i + 1 != data.length) {
            data[i+1][j] = data[i+1][j]/old_value;  // normailze this
          }
        }
        i++; // increment one more time
      }
    }
    else {
      System.err.println("Unkown data type.  Unable to normalize sensor data");
    }
    
    return data;
  }
  
  // Method to get the optical density at wavelength xxx nm
  private float getOpticalDensity(SpectraData abSD, float exWL, float nWL) {
    float od = 1;
    float[][] data = abSD.getData();
    float precision = data[1][0] - data[0][0];
    precision += precision*0.5f;
    
    // find the normalization factor
    float sub_value = 0.0f;
    float nvalue = 0.0f;
    
    if(nWL != 0.0f) { // if normalization wavelength is less than zero then just return intensity
      for(int i = 0; i < data.length; i++) {
        sub_value = Math.abs(nWL - data[i][0]);
        if(sub_value <= precision) {
          nvalue = data[i+1][1];  // possible index out of bounds exception
          break;
        }
      }
    }
    
    // get the normalization value
    for(int i = 0; i < data.length; i++) {
      sub_value = Math.abs(exWL - data[i][0]);
      
      if(sub_value <= precision) {
        od = Math.abs(data[i+1][1] - nvalue); // possible index out of bounds exception
        break;
      }
    }
    
    /*debug code
    System.out.println("Caculating optical density at " + nWL + 
                                    " precision " + precision + " excitation wl " + exWL +
                                    " optical density " + od); //*/
    return od;
  }
  
  // Method to get the area under the curve. It uses the trapezoidal method for integration
  private float getAreaUnderCurve(SpectraData emSD, float minWL, float maxWL) {
    float[][] data = emSD.getData();
    return getAreaUnderCurve(data, minWL, maxWL);
  }
  
  // Method to get the area under the curve. It uses the trapezoidal method for integration
  private float getAreaUnderCurve(float[][] data, float minWL, float maxWL) {
    float area = 0.0f;
    
    float dx =(data[1][0] - data[0][0]);
    
    for( int i=0; i < data.length; i++) {
      float wl = data[i][0];
      if(wl >= minWL && wl <= maxWL && (i + 1) < data.length) {
        float leftSide = data[i][1];
        float rightSide = data[i + 1][1];
        area += dx*(leftSide + rightSide) / 2.0;
      }
    }
    
    return area;
  }
  
  /*
   Method below this point are used to calculate energy transfer value
   */
   
   // Method to calculate the energy transfer efficiency
   public double getEnergyTransferEfficiency(float tau_da, float tau_d) {
     double result = (1 - tau_da/tau_d)*100.0d;
     
     return result;
   }
   
   // Method to calculate the energy transfer probability
   public double getEnergyTransferProbability(float tau_d, float int_d, float int_da) {
     double result = (1/tau_d)*(int_d/(int_da - 1));
     
     return result;
   }
   
   // Method to calculate the Forster distant Ro
   // k = kapp, n = index of refraction, QY = quantum yield, minWL, maxWL,
   // FW = formula weight, C = concentration, L = pathlength
   // matrix = environment of donor and acceptor
   public double getRo(float k2, float n, float QY, float minWL, float maxWL, float FW,
   float C, float L, float[][] em_data, float[][] ab_data, String matrix) {
     double R0 = 0.0d;
     
     /* first get the extinction coefficient array for the absorption data*/
     
     float conc = 0.0f; // concerntration in moles/L
     
     if(matrix.equals("sol-gel")) {
       // must get concerntration for using my sol-gel. Assuming C is in g/ml sol-gel.
       float g_L = C*4.0f*1.0e3f;  // volume decreases by 75% (1/0.25 = 4) in sol-gelP123
       conc = g_L/FW;
     }
     else { // work in solution concerntration is alread correct (moles/L)
       conc = C;
     }
     
     // now convert the absorption values to absorption coefficient or extinction coefficient
     float[][] ec_data = new float[ab_data.length][2];
     for(int i = 0; i < ec_data.length; i++) {
       ec_data[i][0] = ab_data[i][0];
       ec_data[i][1] = ab_data[i][1]/(L*conc);
       //System.out.println("WL : " + ac_data[i][0] + " A : " + ab_data[i][1] + " Epsilon : " + ac_data[i][1]);
     }
     
     /* now calculate the normalize overlap integral J */
     
     // find the index of the minium wavelength index the absorbtion data
     int index = 0;
     
     for(int i = 0; i < ec_data.length; i++) {
       float wl = ec_data[i][0];
       
       if(wl >= minWL) {
         index = i;
         break;
       }
     }
     
     float overlap = 0.0f;
     float norm_overlap = 0.0f;
     
     // caculate the overlap integral now
     float dx = 1.0f;
     
     for(int i = 0; i < (em_data.length -1); i++) {
       float cm1 = em_data[i][0]*1e-7f; // convert wavelength to cm from nm
       float cm2 = em_data[i + 1][0]*1e-7f;
       
       dx = cm2 - cm1;
       
       float x1 = em_data[i][1]*cm1; // multiply by wavelenght to get the quanta
       float x2 = em_data[i + 1][1]*cm2;
       overlap += dx*(x1 + x2)/2.0f;
     }
     
     // calculate the normalized integral
     for(int i = 0; i < (em_data.length -1); i++) {
       if(index < ec_data.length - 1 && em_data[i][0] <= maxWL) {
         float cm1 = em_data[i][0]*1e-7f;
         float cm2 = em_data[i + 1][0]*1e-7f;
         
         dx = cm2 - cm1;
         
         float wl1 = (float)Math.pow(cm1, 4);
         float wl2 = (float)Math.pow(cm2, 4);
         
         float nFd1 = (em_data[i][1]*cm1)/overlap; // normalized emission data
         float nFd2 = (em_data[i + 1][1]*cm2)/overlap;
         
         float x1 = nFd1*ec_data[index][1]*wl1;
         float x2 = nFd2*ec_data[index + 1][1]*wl2;
         
         norm_overlap+= dx*(x1 + x2)/2.0f;
         index++;
       }
     }
     
     /* now calculate the final result */
     float n4 = (float)Math.pow(n, -4.0);
     float f1 = k2*n4*QY*norm_overlap;
     float exp = 1.0f/6.0f;
     float f2 = (float)Math.pow(f1, exp);
     R0 = 9.78e3f*f2;
     
     //* debug code
     System.out.println("DX " + dx);
     System.out.println("QY " + QY);
     System.out.println("N^-4 " + n4);
     System.out.println("K^2 " + k2);
     System.out.println("Overlap " + overlap);
     System.out.println("Normalized Overlap " + norm_overlap);
     System.out.println("f1 " + f1);
     System.out.println("f2 " + f2);
     System.out.println("Exp " + exp);
     System.out.println("R0 in Angs " + R0); //*/
     
     return R0;
   }
   
   // Method to calculate the distance between the donor and the acceptor
   public double getEnergyTransferDistance(float Ro, float E) {
     double value = (1.0d - E)/E;
     double result = Ro*Math.pow(value, 1.0d/6.0d);
     
     return result;
   }
   
   // Method to calculate energy transfer rate
   public double getEnergyTransferRate(float tau_d, float Ro, float r) {
     double result = (1.0f/tau_d)*Math.pow((Ro/r), 6.0d);
     
     return result;
   }
   
}//-- end class EmissionView
