MAMA from MESA website - Wealth Lab Forum

karla2010

#1

9/24/2013 12:32 PM

I coded Ehlers Easy Language version of MAMA from the MESA website, line by line into C#(WLP).

The result matches very closely to the WLP version but differs around rapid changes. The WLP version has smoother edges in the stair steps that one expects from the discontinuous phase transitions in the arctan function.

Is WLP using the same cutoffs as Ehlers for Period and Phase? Is WLP doing some additional smoothing?

Basically all the stairstep locations match within a sample or two but the WLP version has rounded transitions.

Eugene

#2

9/24/2013 1:31 PM

Sometimes indicator authors revise their creations. Not assuming that this is the case here, but MAMA was coded against the the September 2001 issue of Stocks & Commodities magazine. Are you sure the website reflects that version of MAMA?

karla2010

#3

9/24/2013 2:27 PM

ok I don't have a subscription to S&C but I did go to the website and the Traders Tips for 9/2001. One of the contributors showed their code for MAMA and I went through it line by line and it looks the same as the Ehlers version I have.

So the version I have matches the 2001 edition.

Something nonetheless is different. There are some things I like about the WLP version and some things I like about my version.

I'd like to figure out the difference. Either you can release the WLP version so I can go over it and figure out the difference or I'd be happy to publish my version here and you can go over it and decide what is different. I know you are busy and as I said I'd be happy to do the work.

Eugene

#4

9/24/2013 4:07 PM

Sorry, MAMA is a standard indicator. This implies that its implementation is closed source. If you post your script here, we could run it to compare what kind of difference is there.

karla2010

#5

9/24/2013 6:02 PM

ok, here is my version.

The main difference is that I initialize the MAMA computation to get a smooth run in for the values as opposed to the on/off effect seen in WLP.

Ehlers also uses Easy Language so the Atan function is in degrees. I do a conversion from rads to degrees to be consistent with his cutoff(C# uses rads).

I also use double[7] in place of series since there is no need for the detrend, quad and inphase to be series, since only 7 values are needed at any one time and we just dispose of them at the end of the computation.

The remainder follows MESA's published version. I've documented the code the best I could without getting too wordy. The output is of type double[Bars.Count] which I then load into a series for plotting etc.

CODE:
Please log in to see this code.

Size:

Color:

ok, here is my version.

The main difference is that I initialize the MAMA computation to get a smooth run in for the values as opposed to the on/off effect seen in WLP.

Ehlers also uses Easy Language so the Atan function is in degrees. I do a conversion from rads to degrees to be consistent with his cutoff(C# uses rads).

I also use double[7] in place of series since there is no need for the detrend, quad and inphase to be series, since only 7 values are needed at any one time and we just dispose of them at the end of the computation.

The remainder follows MESA's published version. I've documented the code the best I could without getting too wordy. The output is of type double[Bars.Count] which I then load into a series for plotting etc.

[code]
			public double[] cMAMA(DataSeries ds,double ef,double es)
			{
				int ibar,kbar;								//loop indices
				
				double f0,f2,f4,f6;							//HT filter coefficients
				
				DataSeries s0;								//input series

double est,eft;								//input fast, slow alpha
				
				double[] es0 = new double[ds.Count];		//smoothed input
				double[] smm0 = new double[ds.Count];		//MAMA output		
			
				double[] dtrn = new double[7];				//detrender
				double[] sQ1 = new double[7];				//quadrature
				double[] sI1 = new double[7];				//inphase
			
				double sI21,sQ21,sRe1,sIm1,sP1,sPhs1,sMM1;	//unit lag variables
				double sI20,sQ20,sRe0,sIm0,sP0,sPhs0,sMM0;	//current variable
			
				double sJI,sJQ,sdPhs,se0;					//other current variables
			
				double twpi,dc,dc0;							//temp variables common multipliers 
				
				s0 = ds;
				est = es;
				eft = ef;

//initialize filter coefficients
				f0 = 0.0962;	
				f2 = 0.5769;			
				f4 = -f2;
				f6 = -f0;			
			
				//initialize period and HT scale factor
				twpi = 2.0*Math.PI;
				dc0 = 360.0/twpi;
				sP1 = 28.0;
				dc = 0.075*sP1 + 0.54;
			
				//smooth input
			
				for(ibar=3; ibar<ds.Count; ibar=ibar+1)
				{
					es0[ibar] = 0.1*(4.0*s0[ibar]+3.0*s0[ibar-1]
						+2.0*s0[ibar-2]+s0[ibar-3]);
				}
				es0[0] = s0[0];
				es0[1] = 0.33333*(s0[0]+2.0*s0[1]);
				es0[2] = 0.16667*(s0[0]+2.0*s0[1]+3.0*s0[2]);

//initialize detrend
				dtrn[0] = 0.0;
				dtrn[1] = f0*es0[1]+f2*es0[0]+f4*es0[0]+f6*es0[0];
				dtrn[1] = dc*dtrn[1];				
				dtrn[2] = f0*es0[2]+f2*es0[0]+f4*es0[0]+f6*es0[0];
				dtrn[2] = dc*dtrn[2];
				dtrn[3] = f0*es0[3]+f2*es0[1]+f4*es0[0]+f6*es0[0];
				dtrn[3] = dc*dtrn[3];			
				dtrn[4] = f0*es0[4]+f2*es0[2]+f4*es0[0]+f6*es0[0];
				dtrn[4] = dc*dtrn[4];			
				dtrn[5] = f0*es0[5]+f2*es0[3]+f4*es0[1]+f6*es0[0];
				dtrn[5] = dc*dtrn[5];

//initialize quadrature
				sQ1[0] = 0.0;
				sQ1[1] = f0*dtrn[1]+f2*dtrn[0]+f4*dtrn[0]+f6*dtrn[0];
				sQ1[1] = dc*sQ1[1];	
				sQ1[2] = f0*dtrn[2]+f2*dtrn[0]+f4*dtrn[0]+f6*dtrn[0];
				sQ1[2] = dc*sQ1[2];			
				sQ1[3] = f0*dtrn[3]+f2*dtrn[1]+f4*dtrn[0]+f6*dtrn[0];
				sQ1[3] = dc*sQ1[3];			
				sQ1[4] = f0*dtrn[4]+f2*dtrn[2]+f4*dtrn[0]+f6*dtrn[0];
				sQ1[4] = dc*sQ1[4];			
				sQ1[5] = f0*dtrn[5]+f2*dtrn[3]+f4*dtrn[1]+f6*dtrn[0];
				sQ1[5] = dc*sQ1[5];

//initialize inphase
				sI1[0] = 0.0;
				sI1[1] = 0.0;
				sI1[2] = 0.0;
				sI1[3] = dtrn[0];
				sI1[4] = dtrn[1];
				sI1[5] = dtrn[2];
			
				//initialize phase advance for I1 and Q1 by 90 deg
				sJI = f0*sI1[5]+f2*sI1[3]+f4*sI1[1]+f6*sI1[0];
				sJI = dc*sJI;
				
				sJQ = f0*sQ1[5]+f2*sQ1[3]+f4*sQ1[1]+f6*sQ1[0];
				sJQ = dc*sJQ;
				
				//initialize phasor addition for 3 bar averaging
				sI20 = sI1[5] - sJQ;
				sQ20 = sQ1[5] + sJI;
			
				sI21 = sI20;
				sQ21 = sQ20;
			
				//initialize homodyne discriminator
			
				sRe0 = sI20*sI21 + sQ20*sQ21;
				sIm0 = sI20*sQ21 - sQ20*sI21;
			
				sRe1 = sRe0;
				sIm1 = sIm0;
			
				for(ibar=0; ibar<6; ibar=ibar+1)
				{
					smm0[ibar] = es0[ibar];					
				}
			
				//initialize remaining variables
				sP0 = sP1;
				sPhs0 = 0.0;
				sPhs1 = 0.0;
				sMM1 = es0[5];
			
				//start core computation for MAMA
			
				for(ibar=6; ibar<ds.Count; ibar=ibar+1)
				{
					dc = 0.075*sP1 + 0.54;
				
					dtrn[6] = f0*es0[ibar]+f2*es0[ibar-2]+f4*es0[ibar-4]+f6*es0[ibar-6];
					dtrn[6] = dc*dtrn[6];

//inphase and quadrature components
					sQ1[6] = f0*dtrn[6]+f2*dtrn[4]+f4*dtrn[2]+f6*dtrn[0];
					sQ1[6] = dc*sQ1[6];
				
					sI1[6] = dtrn[3];
				
					//advance the phase of I1 and Q1 by 90 deg
					sJI = f0*sI1[6]+f2*sI1[4]+f4*sI1[2]+f6*sI1[0];
					sJI = dc*sJI;
				
					sJQ = f0*sQ1[6]+f2*sQ1[4]+f4*sQ1[2]+f6*sQ1[0];
					sJQ = dc*sJQ;
				
					//phasor addition for 3 bar averaging
					sI20 = sI1[6] - sJQ;
					sQ20 = sQ1[6] + sJI;
				
					//smooth the I and Q components before applying the discriminator
					sI20 = 0.2*sI20 + 0.8*sI21;
					sQ20 = 0.2*sQ20 + 0.8*sQ21;
				
					//homodyne discriminator
					sRe0 = sI20*sI21 + sQ20*sQ21;
					sIm0 = sI20*sQ21 - sQ20*sI21;
				
					sRe0 = 0.2*sRe0 + 0.8*sRe1;
					sIm0 = 0.2*sIm0 + 0.8*sIm1;

//period estimate
					if(sIm0!=0.0 & sRe0!=0.0)
					{
						sP0 = twpi/Math.Atan(sIm0/sRe0);
					}
				
					if(sP0>(1.5*sP1))
					{
						sP0 = 1.5*sP1;					
					}	
				
					if(sP0<(0.67*sP1))
					{
						sP0 = 0.67*sP1;					
					}
				
					if(sP0<6.0)				
					{
						sP0 = 6.0;		
					}	
				
					if(sP0>50.0)
					{
						sP0 = 50.0;	
					}	
				
					sP0 = 0.2*sP0 + 0.8*sP1;
				
					//phase estimate
					if(sI1[6]!=0.0)
					{
						sPhs0 = Math.Atan(sQ1[6]/sI1[6]);
						sPhs0 = dc0*sPhs0;					//convert rads to degs
					}
				
					sdPhs = sPhs1-sPhs0;
				
					if(sdPhs<1.0)
					{
						sdPhs = 1.0;
					}				
				
					//dynamic alpha estimate
					se0 = eft/sdPhs;							//fast alpha
				
					if(se0<est)
					{
						se0 = est;								//slow alpha limit
					}
				
					//MAMA update
					sMM0 = sMM1 + se0*(s0[ibar]-sMM1);
				
					smm0[ibar] = sMM0;
				
					//apply unit delay to all variables for next ibar
				
					for(kbar=0; kbar<6; kbar=kbar+1)
					{
						dtrn[kbar] = dtrn[kbar+1];
						sQ1[kbar] = sQ1[kbar+1];
						sI1[kbar] = sI1[kbar+1];				
					}
				
					sI21 = sI20;
					sQ21 = sQ20;
					sRe1 = sRe0;
					sIm1 = sIm0;
					sP1 = sP0;
					sPhs1 = sPhs0;
					sMM1 = sMM0;
				
				}			
				
				return smm0;
			}
[/code]

Eugene

#6

9/25/2013 1:15 AM

Thanks for sharing your version. Your code looks much longer than Fidelity's implementation in WL6 and Ehlers' EasyLanguage code, too. That's enough to explain the difference.

karla2010

#7

9/25/2013 1:23 AM

hmmmm! My code from the >5 bar is exactly the same as Ehlers, so not sure what you are saying. In any case I'll stick with Ehlers version.

Eugene

#8

9/25/2013 1:55 AM

The Ehlers' version from the PDF above is only 80 lines long. I can estimate Wealth-Lab's version as slightly over 100 lines. Your routine has 230 lines.

karla2010

#9

9/25/2013 2:17 AM

guess we are talking past each other a bit here. If I split a long line, add blank lines, or precompute a piece (e.g. (0.075*period(1) + 0.54) which is used repeatedly, you consider it a complication, I just think of it as efficient and readable code. In any case as I said, for bars > 5, my code is the same as Ehlers core solution. The DSP stuff is exactly the same as Ehlers and I have enough DSP background that I understand what he did, I'm not sure what WLP has done, so I'm just more comfortable with what I can see. I understand the decision to keep certain code proprietary.