#*******************************************************************************************************************
#                       Dispersion Calculator for FIMMWAVE by Klearchos Chaloulos
#                      ---------------------------------------------------
#                                 Copyright Photon Design 2006
#
# This script calculates the dispersion spectrum of the specified node.
# It can do this for one or many modes simultaneously, and for every solver.
# What it does is to calculate the mode list for the waveguide of interest at different wavelengths.
# Then the group index and dispersion can be calculated by applying the appropriate formulas
#
#******************************************************************************************************************



from string import *
from cmath import *
from array import *
from pdPythonLib import *


class matrix:
    def __init__(self, dataType,rows,cols):
        self.m_data=array(dataType,cols*rows*[0])
        self.m_rows=rows
        self.m_cols=cols
    def get(self,row,col):
        a=self.m_data[row*self.m_cols+col]
        return a
    def set(self,row,col,value):
        self.m_data[row*self.m_cols+col]=value
        return
    def rows(self):
        return self.m_rows
    def cols(self):
        return self.m_cols
    def strLine(self,row,AddNewLine):
        strTemp=""
        for j in range(0,self.m_cols):
                strTemp=strTemp+repr(self.m_data[row*self.m_cols+j])+'\t'                
                continue
        if AddNewLine:
            strTemp=strTemp+'\n'
        return strTemp
    def __repr__(self):
        strTemp=""
        for i in range(0,self.m_rows):
            for j in range(0,self.m_cols):
                strTemp=strTemp+repr(self.m_data[i*self.m_cols+j])+'\t'                
                continue
            strTemp=strTemp+'\n'
            continue
        return strTemp
       
def ContainsStr(container,containee):
       a = find(container,containee)
       if a==-1:
           return 0
       return 1

def FindSpaces(str):
       indexArray=array('l')
       length=len(str)
       i=0
       j=0
       jprev=0
       j=find(str," ")
       while j>0 or j==0:
              indexArray.append(j+jprev)
              jprev=j+jprev+1              
              str=str[j+1:]
              j=find(str," ")
              continue
       return indexArray


def InitialiseModeListFMM(startWavelength,nModes):
       #For FMM solvers. Automatically initialise and calculate the mode at starting wavelength
       #Useful before polishing
       global f
       #SolverId=f.Exec("mmainComp.evlist.svp.solvid")
       strBuff=f.Exec("mainComp.evlist.svp.buff")
       #f.Exec("mainComp.evlist.svp.solvid=6")
       nSpaces=FindSpaces(strBuff)
       strBuffNew=strBuff[0:nSpaces[0]+1]+"1"+strBuff[nSpaces[1]:] #starting with 1 1D mode???!!!
       f.Exec("mainComp.evlist.svp.buff="+strBuffNew)
       a=f.Exec("mainComp.evlist.mlp.maxnmodes=1")
       f.Exec("mainComp.evlist.update")
       
       #f.Exec("mainComp.evlist.svp.solvid="+repr(SolverId))
       f.Exec("mainComp.evlist.svp.buff="+strBuff)
       f.Exec("mainComp.evlist.svp.lambda="+repr(startWavelength))       
       a=f.Exec("mainComp.evlist.mlp.maxnmodes="+repr(nModes))
       #print("maxmodes",a)
       #a=f.Exec("mainComp.evlist.unalloc")
       #print("unalloc",a)
       a=f.Exec("mainComp.evlist.update")
       #print ("update",a)
       a=f.Exec("mainComp.evlist.list")
       #print ("list",a)
       a=f.Exec("mainComp.evlist.list[1].neff")
       print("neff",a)
       return

def DoNeffScanPolish():
       #Find the effective index at each step by polishing the mode of the previous step.
       #Suitable for FMM solvers
       global arrWaves
       global f
       global matNeff
       global nSteps
       global nModes
       for i in range(0,nSteps, 1):
             print("Step ",i+1)
             f.Exec("mainComp.evlist.svp.lambda="+repr(arrWaves[i]))
             f.Exec("mainComp.evlist.polishevs")
             for j in range(0,nModes):
                    ctemp=f.Exec("mainComp.evlist.list["+repr(j+1)+"].neff")
                    matNeff.set(i,j,ctemp.real)
                    print ("Mode "+repr(j+1)+", Effective index ="+repr(matNeff.get(i,j)))
                    continue
             print("Done")

             #a=f.Exec("mainComp.evlist.list[1].modedata.tefrac")
             #print ("TE fraction after polish ="+repr(a))
             #ctemp=f.Exec("mainComp.evlist.list[1].neff")
             #arrNeff.append(ctemp.real)
             #print ("Step", i+1, " completed. Effective index =", arrNeff[i])

             continue
       return

def DoNeffScan():
       #Find the effective index at each point by directly finding the mode at each point
       #No polishing. Suitable for FEM and effective index solvers
       global arrWaves
       global f
       global arrNeff
       global nSteps
       global nModes
       a=f.Exec("mainComp.evlist.mlp.maxnmodes="+repr(nModes))
       for i in range(0,nSteps, 1):
             print("Step ",i+1) 
             print("Step ",i+1)
             f.Exec("mainComp.evlist.svp.lambda="+repr(arrWaves[i]))
             f.Exec("mainComp.evlist.update")
             for j in range(0,nModes):
                    ctemp=f.Exec("mainComp.evlist.list["+repr(j+1)+"].neff")
                    matNeff.set(i,j,ctemp.real)
                    print ("Mode "+repr(j+1)+", Effective index ="+repr(matNeff.get(i,j)))
                    continue
             print("Done") 
             continue
       return

def FindNeffSpectrum():
       #Call Fimmwave to evaluate the effctive index at each step. 
       #Call appropriate function to do this according to solver type
       #and polishing options
       global startWavelength
       global nModes
       CurrentSolverType=f.Exec("mainComp.evlist.svp.solvid")
       if (CurrentSolverType == 65 or CurrentSolverType == 73 or
           CurrentSolverType ==17 or CurrentSolverType ==33 or
           CurrentSolverType ==25 or CurrentSolverType ==41):
              print ("Polish of eigenvalues is available with the FMM solver")
              print ("Polishing is much faster, but it could cause problems with many modes")
              print("esp with more than one modes")
              a=raw_input("Polish eigenvalues at each step?(y/n)>>>>")
              if a=="y":
                     print ("Polishing Eigenvalues at each step")                     
                     a=raw_input("Initialise (y/n?)>>>>")
                     if a=="y":
                            InitialiseModeListFMM(startWavelength,nModes)
                     DoNeffScanPolish()
              else:
                     print ("Not polishing, doing full build at eact step")
                     DoNeffScan()              
       else:
              print("Not FMM Type Solver. Will do full build at each step")
              DoNeffScan()
       return      


def FindNg():
       #Find the group index at each step
       global arrWaves
       global arrNeff
       global matNg
       global fWStep
       global nSteps
       global nModes
       for i in range(0,nSteps-1, 1):
             for j in range(0,nModes):
                    der=(matNeff.get(i+1,j)-matNeff.get(i,j))/fWStep
                    matNg.set(i,j,matNeff.get(i,j)-arrWaves[i]*der)
                    continue
             continue
       return

def FindDispersion():
       #Calculate the dispersion spectrum in ps/(nm.Km)
       global arrWaves
       global matNeff
       global fWStep
       global nSteps
       global matDisp
       for i in range(1,nSteps-1, 1):
             for j in range(0,nModes):
                    der2=(matNeff.get(i+1,j)+matNeff.get(i-1,j)-2*matNeff.get(i,j))/pow(fWStep,2)
                    matDisp.set(i-1,j,-arrWaves[i]*10000*der2/3)
                    continue
             continue
       return

def CheckIfFileExists(strFileName):
    #returns a list ret. ret[1] is a string containing a file name. If ret[0]==0, the program can write on this file
    #as a file with that name does not exist or the user decided to overwrite it.
    #If ret[0]==1, the the user decided to abort the procedure
    ret=list()
    ret.append(0)
    ret.append(strFileName)
    strNewFileName=strFileName
    try:
        ff=file(strFileName,'r')
        q=1 #file exists
        ff.close()
    except:
        q=0 #file does not exist
    if q==1:
           print("File already exists. It will be deleted if you proceed.")
           print ("Do you want to specify another filename?(y/n)")
           d=raw_input(">>>>")
           if d=='y':
               strNewFileName=raw_input("Specify full path and filename>>>>")
               r=CheckIfFileExists(strNewFileName)
               q=r[0]
           else:
               q=0
               strNewFileName=strFileName;
    ret[0]=q
    ret[1]=strNewFileName
    return ret
                       

def WriteToFile(strFileName):
       #Write the data (wavelength at each step, effective index, group index and
       #dispersion in a file with name strFileName
       global arrWaves
       global matNeff
       global matNg
       global matDisp
       global nSteps
       global strTemp

       #Check if file already exists. If no exception is raised, file exists.
       r=CheckIfFileExists(strFileName)
       q=r[0];
       strFileName=r[1]
       if q==0:
           print("Writing in file "+strFileName)
           fD=file(strFileName,'w+')
           strTemp=repr(arrWaves[0])+'\t'+matNeff.strLine(0,0)+matNg.strLine(0,0)+'\n'
           fD.write(strTemp)
           fD.flush()
           for i in range(1,nSteps-1, 1):
               strTemp=repr(arrWaves[i])+'\t'+matNeff.strLine(i,0)+matNg.strLine(i,0)+matDisp.strLine(i-1,0)+'\n'
               fD.write(strTemp)
               fD.flush()
               continue
           strTemp=repr(arrWaves[nSteps-1])+'\t'+matNeff.strLine(nSteps-1,0)+'\n'
           fD.write(strTemp)
           fD.flush()
       else:
           print("file operation aborted by user")
                
       return
       

#main program
print ""
print "*****************************************************************"
print "*             Dispersion Calculator for FIMMWAVE                *"
print "*          -----------------------------------                  *"
print "*        by Klearchos Chaloulos (c) Photon Design 2006          *"
print "*                                                               *"
print "*****************************************************************"
print ""

    
f=pdApp()
f.ConnectToApp()
print "Enter the node name of the component."
print "For example if I have a component called comp1"
print "which is stored in the root node under a project called 'project'"
print "then I would enter 'project/comp1':   "
nodeName=raw_input(">>>>")
f.AddCmd("Ref& mainComp=app.findnode({nodeName})")
q = f.Exec("mainComp")
if q==None:
    raise RuntimeError,"Component doesn't exist"
strCmd="mainComp.objtype"
strCompType=f.Exec(strCmd)
if (ContainsStr(strCompType,"RWGuideNode") or ContainsStr(strCompType,"MWGuideNode")
    or ContainsStr(strCompType,"FWGuideNode")):
    print ("Component is "+strCompType+". OK to continue")    
else:
    print ("Component is {strCompType}. Cannot calulate Dispersion for this component.")
    raise RuntimeError,"Component Unsuitable"
strTemp=raw_input("Set Number of modes you are interested in>>>>")
nModes=int(strTemp)
strTemp=raw_input("Set start wavelength>>>>")
startWavelength=float(strTemp)
strTemp=raw_input("Set end wavelength>>>>")
endWavelength=float(strTemp)
strTemp=raw_input("Set number os steps>>>>")
nSteps=int(strTemp,10)
print ("received parameters are:")
print ("Start wavelength:", startWavelength)
print ("End wavelength:", endWavelength)
print ("Number of steps:", nSteps)

fWStep=(endWavelength-startWavelength)/(nSteps-1)
print ("Wavelength Step=", fWStep)


arrWaves=array('d')
matNeff=matrix('d',nSteps,nModes)
matNg=matrix('d',nSteps-1,nModes)
matDisp=matrix('d',nSteps-2,nModes)

for i in range(0,nSteps):
       fTemp=startWavelength+i*fWStep
       arrWaves.append(fTemp)

#Find Neffs
FindNeffSpectrum()
print (matNeff)

FindNg()
print("Group index:")
print(matNg)

FindDispersion()
print("Dispersion:")
print(matDisp)

a=raw_input("Write to file?(y/n)>>>>")
if a=="y":
       strFileName=raw_input("Specify full path and filename>>>>")
       WriteToFile(strFileName)
      


#delete variables
del nodeName
del arrWaves
del matNeff
del matDisp
del matNg
del strTemp
del strCmd
del strCompType
del q
del i
del startWavelength
del endWavelength
del nSteps
del fWStep
del fTemp
del f