import java.awt.Color;
import java.awt.Dimension;
import java.awt.Image;
abstract public class MultiPathConvolution extends Content{
public MultiPathConvolution(Dimension dim, Image img, String problemName){
super(dim,img,problemName);
}
////////////////////////////////////////////////////
// AFFECTS CHILD CLASSES ONLY WHEN NOT OVERRIDDEN.
//==================================================
protected void spawnParametersAfterUserApplied(){
spanParsInMPLevel();
displySpawnedParsInMPLevel();
simulationBlocked=false;
simulateThresholding();
}
////////////////////////////////////////////////////
//=============================
// Sub Model
//-----------------------------
int subModel; //Index of following subModels:
String[] subModelTitle={ "OOK and PAM", "DAPPM", "PPM", "PAPM"};
int shotNoisePresented=1;
int amInteferenceSummationPoints=2;
//=============================
//=============================
// Default parameters
// See description at line:
// User Input Prompts:
//-----------------------------
double ceilingHeight=3.5;
int Amax=3;
int L=1;
double Rb=100000000;
double SNR=7;
double amSAR=1.0;
double amInterferencePeriodTi=25.0E-6;
//=============================
//=============================
//Derivative parameters:
// See description at line:
// To display this strings in graph background:
//-----------------------------
double a;
//SNR not in dB form:
double SN;
double T;
double avLength;
double aphabetCount;
double M;
double bitsPerChip;
double scaled_chip_length;
int tapsNumber;
double[] beta; //Discretized h.
//min non-zero Intensity/average Intensity:
double lambda;
//Submodel specific:
double OOK_threshold=0.5;
//=============================
//Simulation control:
protected boolean simulationBlocked;
//Intermediate Results:
//double B; //BER
double BAmb=0.0; //BER over Ti - Ambient interference period.
/////////////////////////////////////////////////////////////////
// Assign Parameters
//================================================================
//================================================================
//Can be used in child classes directly.
//If not used in child classes, then do not affect them.
//----------------------------------------------------------------
final protected void initPreconstructorInMPLevel(){
modelDescription=new String[200];
modelDescription[0]="Program: Stub for MultiPath Indoor Dispersion.";
tapsSimulationLimit=31;
amPrepare();
dmnStartF=0.0;
dmnRangeF=8.0;
dmnRangeFDown=1.0;
dmnStartX=0.0;
dmnRangeX=3.0;
grPoints=400;
grStartF=600;
functionCOUNT=7;
functionColor=new Color[]{
new Color(255,0,0),
new Color(220,100,0),
new Color(0,0,200),
new Color(150,0,150),
new Color(0,150,0),
new Color(0,255,0),
new Color(100,100,100,50),
new Color(200,0,0,50),
new Color(222,0,200,50)
};
functionTitle=new String[]{
"h",
"beta",
"Origin",
"Loss",
"Origin",
"Gain",
"FirstLossThreshold",
"Ambient V",
"Ambient V over chip"
};
}
//To be explicitly used in child classes:
final protected void spanParsInMPLevel(){
//con("Spawning parameters InMPLevel...");
a=2.*ceilingHeight/300000000.0;
setSubModel();
bitsPerChip=M/avLength;
T=bitsPerChip/Rb;
scaled_chip_length=T/a;
//-------------------------------------------------------
///stimation of size of sequence beta:
//- - - - - - - - - - - - - - - - - - - - - - - - - - - -
double accuracyEps=1.0e-3;
double hThresholdTs= Math.exp( Math.log(1.0/accuracyEps)*(1.0/6.0) ) - 1;
if(hThresholdTs<1.0)hThresholdTs=1.0;
double wTapsNumber=hThresholdTs/scaled_chip_length;
tapsNumber = ((int)Math.floor(wTapsNumber)) + 1; //1 is taken for safety.
double wTruncationParameterForCuriosity=wTapsNumber / Math.floor(a/T);
con("wTruncationParameterForCuriosity"+wTruncationParameterForCuriosity);
//Was:
//tapsNumber= ((int)Math.floor(a/T*1.25)) + 1;
//-------------------------------------------------------
//Convert SN from dB to numbers:
SN=Math.exp( SNR/10.0*Math.log(10.0) );
setSequenceMemoryInMPLevel();
assignBetaKToArray();
}
//To display this strings in graph background:
final protected void displySpawnedParsInMPLevel(){
//con("subModelTile[0]="+subModelTitle[0]+" " + subModel);
modelDescription[1]="Sub Model = "+subModelTitle[subModel];
modelDescription[4]="A="+Amax;
modelDescription[5]="L="+L;
modelDescription[6]="SNR="+SNR+",db SNR="+SN;
modelDescription[7]="Symbols in Alphabet="+aphabetCount;
modelDescription[8]="M="+M;
modelDescription[9]="BitRate, Rb="+Rb+", bits/second.";
modelDescription[10]="Average Length="+avLength;
modelDescription[11]="Bits per chip="+bitsPerChip;
modelDescription[12]="Chip Length="+T+", seconds.";
modelDescription[13]="Room Height="+ceilingHeight+", meters.";
modelDescription[14]="Dispersion Length Scale, a="+a+", seconds.";
modelDescription[15]="Lambda="+lambda;
modelDescription[16]="tapsNumber="+tapsNumber;
}
//----------------------------------------------------------------
//Can be used in child classes directly.
//================================================================
//================================================================
//Can not be used in child classes directly.
//----------------------------------------------------------------
private void setSequenceMemoryInMPLevel(){
//con(" Memory");
int len=tapsNumber+L+1;
beta=new double[tapsNumber];
b=new int[tapsNumber+L];
result=new double[len];
//----------------------------------------------------------
//Set memory for signals which demonstrate ISI loss or gain:
inputChipLost=new int[len];
inputChipGained=new int[len];
resultChipLost=new double[len];
resultChipGained=new double[len];
for(int i=0; i<result.length; i++){
inputChipLost[i]=0;
inputChipGained[i]=0;
}
//----------------------------------------------------------
}
private void setSubModel(){
int i=L;
switch(subModel){
case 0: //PAM
avLength=L;
aphabetCount=Amax+1;
while(--i>0){
aphabetCount*=Amax+1;
}
M=Math.log(aphabetCount)/Math.log(2.0);
lambda=2.0/Amax;
break;
case 1: avLength=(1+L)/2; //In Chips.
aphabetCount=Amax*(1+L)*L/2;
M=Math.log(aphabetCount)/Math.log(2.0);
lambda=1.0/((Amax+1.0)/2/avLength);
break;
case 2: //PPM
avLength=L; //In Chips.
aphabetCount=L;
M=Math.log(aphabetCount)/Math.log(2.0);
lambda=L;
break;
case 3: //PAPM
avLength=L; //In Chips.
aphabetCount=L*Amax;
M=Math.log(aphabetCount)/Math.log(2.0);
lambda=2.0*L/(Amax+1);
con("Amax="+Amax+" L="+L+" lambda="+lambda);
break;
}
}
//----------------------------------------------------------------
//Can not be used in child classes directly.
//================================================================
//Assign Parameters
//////////////////////////////////////////////////////////////////
//==========================================
// Algorithms
//------------------------------------------
protected double hh(double q){
double power=(q+1);
double p=power*power;
double p4=p*p;
p*=p4; //6
p*=power; //7
return 1/p*6.0;
}
//Returns area under h-function over interval [k*step,k*step+step]:
//Input: step = scaled chip length.
private double betaPortion(double k, double step){
double t=k*step;
double power=t+1;
double p=power*power;
double p4=p*p;
double value1=1/(p*p4); //6
t+=step;
power=t+1;
p=power*power;
p4=p*p;
return value1-1/(p*p4);
}
private void assignBetaKToArray(){
for(int k=0; k<tapsNumber; k++){
beta[k]=betaPortion((double)k, scaled_chip_length);
//con("beta="+beta[k]);
}
}
//Calculate all result chips:
protected void convolve(int[] inputPul,double[] outputPul, boolean drawNumbers){
convolve(0, inputPul,outputPul, drawNumbers);
}
//Calculate only from start result chip:
protected void convolve(int start, int[] inputPul,double[] outputPul, boolean drawNumbers){
int nOut=outputPul.length;
int nIn=inputPul.length;
for(int i=start; i<nOut; i++){
double s=0;
for(int j=0; j<tapsNumber; j++){
int tail=i-j;
if(tail<0 || tail>=nIn)break;
s+=inputPul[tail]*beta[j];
}
outputPul[i]=s;
if(drawNumbers) {
int drawIndex=25+i;
if(drawIndex<modelDescription.length){
String ss=" conv "+i+" "+s+" beta="+beta[Math.min(beta.length-1,i)];
if(i<nIn)ss+=" input="+inputPul[i];
modelDescription[drawIndex]=ss;
}
}
}
}
//------------------------------------------
// Algorithms
//==========================================
//==========================================
// Simulation
//------------------------------------------
protected int[] b;
//do1 replace with bh:
protected double[] result;
protected int tapsSimulationLimit;
protected void simulateThresholding(){
//Demo:
//Clear display:
int i;
for(i=20; i<modelDescription.length; i++){
modelDescription[i]=null;
}
if(simulationBlocked)return;
if(subModel>0){
modelDescription[20]="BER for sub Model " + subModelTitle[subModel]+" is not implemented.";
simulationBlocked=true;
return;
}
if(tapsSimulationLimit<=tapsNumber){
modelDescription[20]="BER for sub Model " + subModelTitle[subModel]+". Taps Number exceeded limit.";
//Set B to 1.0 to attract attention:
BAmb=1.0;
simulationBlocked=true;
return;
}
//PAM and OOK:
int unitEventsCount=2;
int eventsCount=Amax+1;
int tL=Math.min(tapsSimulationLimit,tapsNumber);
i=tL;
while(--i>0){
unitEventsCount*=2;
eventsCount*=Amax+1;
};
modelDescription[22]="Chip sequences count " + eventsCount +
" Unit sequences count=" + unitEventsCount;
int EVENTS_MEASURE_LIMIT=1000000;
if(eventsCount>EVENTS_MEASURE_LIMIT){
modelDescription[20]="Chip sequences count " + eventsCount +
" > " + EVENTS_MEASURE_LIMIT + " which seems too long to process.";
return;
}
BAmb=0.0;
double amInterferenceStep=amInterferencePeriodTi/amInteferenceSummationPoints;
for(int iXAm=0; iXAm<amInteferenceSummationPoints; iXAm++){
double t=amInterferenceStep*iXAm;
double BB=0.0;
double Bup=0;
double Bdown=0;
double BupISI=0;
double BdownISI=0;
firstLossThreshold=-1;
double firstGainThreshold=-1;
//Shortcuts:
double teta=lambda*OOK_threshold;
double amK=1.0/amSAR;
//con("teta="+teta);
boolean HeavisideMode=0==shotNoisePresented;
int k=tL-1; //Probe slot to count statistics.
for(int e=0; e<unitEventsCount; e++){
//-------------------------------
//Generate signal of units:
int mask=e;
int weight=1;
int count=0;
for(int slot=0; slot<tL; slot++){
int ampl=mask&1;
if(ampl>0){ weight*=Amax; count++; }
b[k-slot]=ampl;
mask>>>=1;
}
//-------------------------------
for(int iW=0; iW<weight; iW++){
//----------------------------------------------
//Decompose iW and assign amplitudes to a signal
//of units:
i=tL;
int weightS=iW;
while(--i>-1){
if(b[i]>0){
int reminder=weightS%Amax;
b[i]=reminder+1;
weightS=(weightS-reminder)/Amax;
}
}
//----------------------------------------------
convolve(k,b,result,false);
//con("e="+e);
//con("iW="+iW);
//con("bh="+result[k]);
double S=lambda*result[k];
//con("S="+S);
double Z=S; //no ambient light yet.
if(amInteferenceSummationPoints>1){
Z=Z+amK*am_vi(t, scaled_chip_length*a);
}
int bk=b[k];
double nUp=lambda*bk+teta-Z; //noise Up
double nDown=Z-(lambda*bk-teta); //noise Down
double up=0;
double down=0;
//We have three principal cases, bk=0, bk=A, bk in the middle.
if(0==bk){
up=UFun.erfh(SN*nUp,HeavisideMode);
}else if(Amax==bk){
down=UFun.erfh(SN*nDown,HeavisideMode);
}else{
up=UFun.erfh(SN*nUp,HeavisideMode);
down=UFun.erfh(SN*nDown,HeavisideMode);
}
if(0==subModel){
if(S>lambda*bk+teta){
if(0.0==BupISI){
convolve(b,result,false);
copyArray(b,inputChipGained);
copyArray(result,resultChipGained);
firstGainThreshold=bk+teta/lambda;
}
BupISI++;
}
if(S<lambda*bk-teta ){
if(0.0==BdownISI){
convolve(b,result,true);
copyArray(b,inputChipLost);
copyArray(result,resultChipLost);
firstLossThreshold=bk-teta/lambda;
}
BdownISI++;
}
}
Bup+=up;
Bdown+=down;
BB+=up;
BB+=down;
//con(" B="+BB+" nUp="+nUp+" up="+up+" Bup="+Bup+" nDown="+nDown+" down="+down+" Bdown="+Bdown);
} //for(int iW=0; iW<weight; iW++){
}// e
//con("events Count="+eventsCount);
BB/=eventsCount;
Bup/=eventsCount;
Bdown/=eventsCount;
BdownISI/=eventsCount;
BupISI/=eventsCount;
//con(" B="+BB+" Bup="+Bup+" Bdown="+Bdown);
//-------------------------------------
//Output:
modelDescription[21]="";
modelDescription[22]="";
if(1==amInteferenceSummationPoints){
modelDescription[20]="BER, BERup, BERdown = "+BB+", "+Bup+", "+Bdown;
if(0==subModel){
modelDescription[21]="BERupISI, BERdownISI = "+BupISI+", "+BdownISI;
modelDescription[22]="firstLossThreshold, firstGainThresho = "+
firstLossThreshold+", "+firstGainThreshold;
}
}
//-------------------------------------
BAmb=BAmb+BB/amInteferenceSummationPoints;
} // for(int iXAm=0; iXAm<amInteferenceSummationPoints; iXAm++){
if(1<amInteferenceSummationPoints)modelDescription[20]="BER "+BAmb;
}
//------------------------------------------
// Simulation
//==========================================
//==========================================
// Draw Graphs
//------------------------------------------
protected double drawSignal(double t, double[] stepSignal){
int k=(int)(t/scaled_chip_length);
if(k>=stepSignal.length)return 0;
return (double) stepSignal[k];
}
//Integer version of the method:
protected double drawSignal(double t, int[] stepSignal){
int k=(int)(t/scaled_chip_length);
if(k>=stepSignal.length)return 0;
return (double) stepSignal[k];
}
abstract protected double functionSwitch(int fIx, double x);
//------------------------------------------
// Draw Graphs
//==========================================
//==========================================
// Ambient light
//------------------------------------------
//%
//%Credit: Data has been used from the following source:
//% Adriano J.C. Moreira, Rui T. Valadas and A.M. de Oliveira Duarte.
//% Optical interference produced by artificial light.
//% Wireless Networks 3 (1997) 131–140
//%
protected double[] am_b;
protected double[] am_c;
protected double[] am_zeta;
protected double[] am_fi;
protected double[] am_d;
protected double[] am_teta;
protected double am_A1_reciprocal;
protected double am_A2_reciprocal;
//Fundamental frequency of high frequency component in (7), Hz:
protected double am_fh;
static final double PI2=2*Math.PI;
protected void amPrepare(){
am_b=new double[21];
am_c=new double[21];
am_zeta=new double[21];
am_fi=new double[21];
am_d=new double[12];
am_teta=new double[12];
am_A1_reciprocal=1.0/5.9;
am_A2_reciprocal=1.0/2.1;
//Fundamental frequency of high frequency component in (7), Hz:
am_fh=37.5E3;
double log10 = Math.log(10);
for(int i=1; i<21; i++){
am_b[i]=Math.exp( log10*( -13.1*Math.log(100*i-50) +27.1 )/20 );
am_c[i]=Math.exp( log10*( -20.8*Math.log(100*i) + 92.4 )/20 );
}
//[7] Table I:
double[] amAux=new double[]{
1, 4.65, 0.00, 11, 1.26, 6.00,
2, 2.86, 0.08, 12, 1.29, 6.17,
3, 5.43, 6.00, 13, 1.28, 5.69,
4, 3.90, 5.31, 14, 0.63, 5.37,
5, 2.00, 2.27, 15, 6.06, 4.00,
6, 5.98, 5.70, 16, 5.49, 3.69,
7, 2.38, 2.07, 17, 4.45, 1.86,
8, 4.35, 3.44, 18, 3.24, 1.38,
9, 5.87, 5.01, 19, 2.07, 5.91,
10, 0.70, 6.01, 20, 0.87, 4.88
};
for(int ii=0; ii<10; ii++){
for(int j=0; j<2; j++){
int pos=ii*6+j*3;
int i=ii+10*j+1;
am_zeta[i]=amAux[pos+1];
am_fi[i]=amAux[pos+2];
}
}
//Check our selves:
for(int i=1; i<=10; i++){
int j=i+10;
//con(""+i+" "+am_zeta[i]+" "+am_fi[i]+" "+j+" "+am_zeta[j]+" "+am_fi[j]);
}
//[7] Table II. Amplitude and phase parameters for high-frequency components.
double[] amAux2=new double[]{
0, -22.22, 5.09, 6, -39.30, 3.55,
1, 0.00, 0.00, 7, -42.70, 4.15,
2, -11.50, 2.37, 8, -46.40, 1.64,
3, -30.00, 5.86, 9, -48.10, 4.51,
4, -33.90, 2.04, 10, -53.10, 3.55,
5, -35.30, 2.75, 11, -54.90, 1.78
};
for(int ii=0; ii<6; ii++){
for(int j=0; j<2; j++){
int pos=ii*6+j*3;
int i=ii+6*j;
am_d[i]=amAux2[pos+1];
am_teta[i]=amAux2[pos+2];
}
}
//Check our selves:
for(int i=0; i<=5; i++){
int j=i+6;
//con(""+i+" "+am_d[i]+" "+am_teta[i]+" "+j+" "+am_d[j]+" "+am_teta[j]);
}
}
//Here we are programming formula [7 (7)]:
//Factors R and Pm are not included into m(t) as in formula (7):
protected double am_mi(double t){
//Calculate member 1 in (7), RPm.
double sum1=1.0;
//Calculate member 2:
double sum2=0.0;
for(int i=1; i<=20; i++){
sum2 += am_b[i]*Math.cos(PI2*(100*i-50)*t+am_zeta[i])+
am_c[i]*Math.cos(PI2*(100*i*t+am_fi[i]));
}
sum2*=am_A1_reciprocal;
//Calculate member 3:
double sum3=am_d[0]*Math.cos(PI2*am_fh*t+am_teta[0]);
for(int i=1; i<=11; i++){
sum3 += am_d[i]*Math.cos(PI2*2*i*am_fh*t+am_teta[i]);
}
sum3*=am_A2_reciprocal;
double sum=sum1+sum2+sum3;
//UTil.con("amV="+sum);
return sum;
}
//Calculate ambient light average energy contribution over time t2-t1 starting
//from time t1.
//Returns: 1/dt*Integral_over amV by dt.
//Purpose: calculate noise applied to one chip:
protected double am_vi(double t1, double dt){
double a1=0;
double delta=0;
double alpha=0;
//Calculate member 1 in (7), RPm.
double sum1=1.0*dt;
//Calculate member 2:
double sum2=0.0;
for(int i=1; i<=20; i++){
//Calculate member 2:
double a0=PI2*100*i;
a1=a0-PI2*50;
double delta0=a0*dt*0.5;
double delta1=a1*dt*0.5;
double alpha1=a1*t1+am_zeta[i]+delta1;
double alpha0=a0*t1+am_fi[i]+delta0;
sum2 = sum2+am_b[i]*2.0*(
Math.sin(delta1)*Math.cos(alpha1)/a1+
Math.sin(delta0)*Math.cos(alpha0)/a0
);
}
sum2=sum2*am_A1_reciprocal;
//Calculate member 3:
a1=PI2*am_fh;
delta=a1*dt*0.5;
alpha=a1*t1+am_teta[0]+delta;
double sum3=am_d[0]*2.0*Math.sin(delta)*Math.cos(alpha)/a1;
for(int i=1; i<=11; i++){
a1=PI2*2*i*am_fh;
alpha=a1*t1+am_teta[i]+delta;
sum3 += am_d[i]*Math.sin(delta)*Math.cos(alpha)/a1;
}
sum3*=am_A2_reciprocal;
double sum=(sum1+sum2+sum3)/dt;
//UTil.con("amV="+sum);
return sum;
}
//------------------------------------------
// Ambient light
//==========================================
//==========================================
// Auxiliary procedures
//------------------------------------------
protected void copyArray(int[] f, int[] t){
int lim=Math.min(f.length,t.length);
for(int i=0; i<lim; i++){
t[i]=f[i];
}
}
protected void copyArray(double[] f, double[] t){
int lim=Math.min(f.length,t.length);
for(int i=0; i<lim; i++){
t[i]=f[i];
}
}
//------------------------------------------
// Auxiliary procedures
//==========================================
//==========================================
// Demonstration parameters.
// To store signals which show failure
// detection of a chip due ISI:
// Used: in simulateThreshold() method.
//------------------------------------------
protected double firstLossThreshold;
protected int[] inputChipLost;
protected int[] inputChipGained;
protected double[] resultChipLost;
protected double[] resultChipGained;
//==========================================
}
Copyright (C) 2009 Konstantin Kirillov