package ppmrmn; import java.awt.*; import java.applet.*; import java.awt.event.*; //1 février 2004 public class OnePowderJDK118 extends Graph implements Runnable, Population { private int numVars = 3, //number of parameters to be fitted with simplex nbou = 1, //1 + number of pulse duration increment ligne = 2, //transition position in the density matrix nbBarre = 25;//nomber of columns in the powder pattern private char choix = 'C'; private boolean stopFit = false, showCurve = true; private double coeff = 5, //for signal normalisation use pop = 2, //for signal normalisation use dureePremiere = 0, //inital pulse duration dureeFin = 1, //final pulse duration dureePas = 1, //increment of the pulse duration amplitude[] = {0.0}, //FID intensity listExpIntensity[] = {0.0}, //experimental intensity /*matEstimation[][] is the matrix for the initial guess values of fitting parameters*/ matEstimation[][] = new double[numVars + 1][numVars]; //°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°° // Dimensions 660LX390H of the applet to be introduced in its WEB page //°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°°° private int hautLB = 20, //label height hautBouton = 30, //button height hautPanel = 390, //panel height population[] = new int[nbBarre]; //for powder pattern //PanelSpin includes panelSpinLabel, panelSpinData, and panelSpinUnit Panel panelSpin = new Panel(); Panel panelSpinLabel = new Panel(); Panel panelSpinData = new Panel(); Panel panelSpinUnit = new Panel(); //panelTransition includes panelSimulation and scrollPaneSimulation Panel panelSimulation = new Panel(); Panel panelTransition = new Panel(); ScrollPane scrollPaneSimulation = new ScrollPane(); //PanelFit includes panelFitEstim, panelFitBouton, and panelBlanc Panel panelFit = new Panel(); Panel panelFitEstim = new Panel(); Panel panelFitBouton = new Panel(); Panel panelBlanc = new Panel(); Choice choiceSpin = new Choice(); Choice choiceEta = new Choice(); Choice choiceTransition = new Choice(); Choice choiceSimplex = new Choice(); Button buttonRun = new Button("Run"); Button buttonFit = new Button("Fit"); Button buttonGraph = new Button("DeactivateGraph"); Button buttonPrevious = new Button("Previous"); Label labelPulse = new Label("Pulse"); Label labelDuree = new Label("Length (µs)"); Label labelLine = new Label("Line"); Label labelIntensite = new Label("Intensity"); Label labelOf = new Label("of"); Label labelExperimental = new Label("Experimental"); Label labelExpIntensite = new Label("Intensity"); Label labelAdvice = new Label("Deactivate graph speeds up fitting!!!"); TextField textFieldOmegaQ = new TextField("0.0"); TextField textFieldOmegaRF = new TextField("50.0"); TextField textFieldDerniereDuree = new TextField("10"); TextField textFieldPas = new TextField("1"); TextField textFieldPremiereDuree = new TextField("0"); TextField textFieldOmegaQEst = new TextField("5"); TextField textFieldOmegaRFEst = new TextField("55"); //100 rows, 7 columns, 2 = SCROLLBARS_NONE TextArea textAreaDuree = new TextArea("", 100, 7, TextArea.SCROLLBARS_NONE); TextArea textAreaIntensite = new TextArea("", 100, 7, TextArea.SCROLLBARS_NONE); TextArea textAreaExpIntensite = new TextArea(12, 7); //experimental intensity TextArea textAreaParameter = new TextArea(5, 7); //fitting parameters /**Construct the applet*/ public OnePowderJDK118() { } /**Initialize the applet*/ public void init() { try { jbInit(); } catch(Exception e) { e.printStackTrace(); } } /**Component initialization*/ private void jbInit() throws Exception { setBackground(Color.white); setLayout(null); buttonGraph.setBounds(75, 355, 150, hautBouton); add(buttonGraph); buttonGraph.setVisible(false); labelAdvice.setBounds(15, 320, 220, hautBouton); add(labelAdvice); labelAdvice.setFont(new Font("TimesRoman", Font.PLAIN, 10)); //necessary for WEB browsers other than IE5.5: labelAdvice.setBackground(Color.pink); labelAdvice.setVisible(false); panelSpin.setBounds(0, 0, 220, hautPanel); add(panelSpin); panelSpinLabel.setBounds(0, 0, 90, hautPanel-3); panelSpin.add(panelSpinLabel); panelSpinData.setBounds(90, 0, 80, hautPanel-3); panelSpin.add(panelSpinData); panelSpinUnit.setBounds(180, 0, 50, hautPanel-3); panelSpin.add(panelSpinUnit); panelTransition.setBounds(220, 0, 220, hautPanel); add(panelTransition); panelFit.setBounds(440, 0, 220, hautPanel); add(panelFit); panelFit.setVisible(false); //**************************** panelSpin ********************************* panelSpin.setBackground(Color.pink); panelSpin.setLayout(null); panelSpin.setFont(new Font("TimesRoman", Font.PLAIN, 12)); choiceSpin.add("3/2"); choiceSpin.add("5/2"); choiceSpin.add("7/2"); choiceSpin.add("9/2"); choiceSpin.select(0); choiceEta.add("0.0"); choiceEta.add("0.1"); choiceEta.add("0.2"); choiceEta.add("0.3"); choiceEta.add("0.4"); choiceEta.add("0.5"); choiceEta.add("0.6"); choiceEta.add("0.7"); choiceEta.add("0.8"); choiceEta.add("0.9"); choiceEta.add("1.0"); choiceEta.select(0); population = POPULATION_ETA00; //8 lines, 1 column, 0 horizontal space, and 10 vertical spaces panelSpinLabel.setLayout(new GridLayout(8, 1, 0, 10)); panelSpinData.setLayout(new GridLayout(8, 1, 0, 10)); panelSpinUnit.setLayout(new GridLayout(8, 1, 0, 10)); panelSpinLabel.setBackground(Color.pink); panelSpinData.setBackground(Color.pink); panelSpinUnit.setBackground(Color.pink); panelSpinLabel.add(new Label("01-02-2004")); panelSpinData.add(new Label("")); panelSpinUnit.add(new Label("")); panelSpinLabel.add(new Label(" Spin")); panelSpinData.add(choiceSpin); panelSpinUnit.add(new Label("")); panelSpinLabel.add(new Label(" RF field")); panelSpinData.add(textFieldOmegaRF); panelSpinUnit.add(new Label("kHz")); panelSpinLabel.add(new Label(" QCC")); panelSpinData.add(textFieldOmegaQ); panelSpinUnit.add(new Label("kHz")); panelSpinLabel.add(new Label(" Eta")); panelSpinData.add(choiceEta); panelSpinUnit.add(new Label("")); panelSpinLabel.add(new Label(" MinLength")); panelSpinData.add(textFieldPremiereDuree); panelSpinUnit.add(new Label("µs")); panelSpinLabel.add(new Label(" MaxLength")); panelSpinData.add(textFieldDerniereDuree); panelSpinUnit.add(new Label("µs")); panelSpinLabel.add(new Label(" Step")); panelSpinData.add(textFieldPas); panelSpinUnit.add(new Label("µs")); //************************** panelTransition ***************************** panelTransition.setBackground(Color.pink); panelTransition.setLayout(null); panelTransition.setFont(new Font("TimesRoman", Font.PLAIN, 12)); labelPulse.setBounds(15, 7, 85, hautLB); panelTransition.add(labelPulse); labelDuree.setBounds(15, 27, 85, hautLB); panelTransition.add(labelDuree); labelLine.setBounds(100, 7, 85, hautLB); panelTransition.add(labelLine); labelIntensite.setBounds(100, 27, 100, hautLB); panelTransition.add(labelIntensite); panelSimulation.setLayout(new GridLayout(1, 2)); panelSimulation.add(textAreaDuree); panelSimulation.add(textAreaIntensite); textAreaDuree.setEditable(false); textAreaIntensite.setEditable(false); scrollPaneSimulation.setBounds(10, 50, 180, 250); scrollPaneSimulation.add(panelSimulation); panelTransition.add(scrollPaneSimulation); labelOf.setBounds(90, 305, 85, hautLB); panelTransition.add(labelOf); choiceTransition.add("Central Transition"); choiceTransition.add("Satellite Transitions"); choiceTransition.add("All the Transitions"); choiceTransition.select(0); choiceTransition.setBounds(10, 325, 180, 30); panelTransition.add(choiceTransition); buttonRun.setBounds(50, 355, 100, hautBouton); panelTransition.add(buttonRun); //**************************** panelFit ********************************** panelFit.setBackground(Color.orange); panelFit.setLayout(null); panelFit.setFont(new Font("TimesRoman", Font.PLAIN, 12)); labelExperimental.setBounds(10, 7, 100, hautLB); panelFit.add(labelExperimental); labelExpIntensite.setBounds(10, 27, 100, hautLB); panelFit.add(labelExpIntensite); textAreaExpIntensite.setBounds(10, 50, 80, 160); panelFit.add(textAreaExpIntensite); textAreaParameter.setBounds(10, 212, 200, 140); panelFit.add(textAreaParameter); panelFitBouton.setBounds(10, 355, 200, hautBouton); panelFit.add(panelFitBouton); panelFitBouton.setLayout(new GridLayout(1, 2)); panelFitBouton.add(buttonPrevious); panelFitBouton.add(buttonFit); buttonPrevious.setEnabled(false); panelFitEstim.setBounds(110, 5, 100, 206); //5 lines, 2 columns, zero horizontal and vertical spaces panelFitEstim.setLayout(new GridLayout(5, 2, 0, 0)); panelFit.add(panelFitEstim); panelFitEstim.setBackground(Color.orange); panelFitEstim.setFont(new Font("TimesRoman", Font.PLAIN, 12)); choiceSimplex.add("1.1"); choiceSimplex.add("1.2"); choiceSimplex.add("1.3"); choiceSimplex.add("1.4"); choiceSimplex.select(2); panelFitEstim.add(new Label("Simplex:")); panelFitEstim.add(choiceSimplex); panelFitEstim.add(new Label("RF field:")); panelFitEstim.add(new Label("")); panelFitEstim.add(textFieldOmegaRFEst); panelFitEstim.add(new Label(" kHz")); panelFitEstim.add(new Label("QCC:")); panelFitEstim.add(new Label("")); panelFitEstim.add(textFieldOmegaQEst); panelFitEstim.add(new Label(" kHz")); //valeurs exemples String s5 = 0+"\n"+1119.8+"\n"+2130.1+"\n"+2931.9+"\n"+3446.6+"\n"+3624+"\n"; s5 += 3446.6+"\n"+2931.9+"\n"+2130.1+"\n"+1119.8+"\n"+0+"\n"; textAreaExpIntensite.setText(s5); }//end of jbInit() private void calculIntensite(double omegaQTmp, double omegaRFTmp) { RealMat matH, matHD, matR, matRt, matIz, matR1, matR2; ComplexMat matHH, matHHC, matC1, matC2, matC3, matC4; //omegaQTmp <-> QCC, and omegaRFTmp in kHz unit Nmr myNmr = new Nmr(omegaQTmp, omegaRFTmp); //omegaQ <-> QCC, and omegaRF in radian/second unit double deltaY = 2.0/nbBarre, //width of a column in the powder pattern spin1 = myNmr.getSpin(), //get the spin value rap = 3.0/(8*spin1*(2*spin1 - 1)), //factor associated with QCC tmp = myNmr.getOmegaQ(); //QCC in radian/second unit for (int i = 0; i <= nbou; i++) amplitude[i] = 0.0; //clear //k < population.length for (int k = 0 ; k < population.length; k++) { myNmr.setOmegaQ(k*deltaY*rap*tmp); //set the value of omegaQ matH = myNmr.dataIn(); //find eigenvectors and eigenvalues of matH matHD = myNmr.jacRot(matH); //matHD is the eigenvalue matrix of matH matR = myNmr.getMatR(); //matR is the eigenvector matrix of matH matRt = RealMat.transpose(matR); //matRt is the transposed matrix of matR matIz = myNmr.matriceIz(); //matrix of spin operator Iz matR1 = myNmr.mult(matIz, matR); //matrix of Iz expressed in the PAS matR2 = myNmr.mult(matRt, matR1); //of matH //action of the pulse for (int j = 0; j <= nbou; j++) { double tp = dureePremiere + j*dureePas; matHH = myNmr.matriceHamilton(matHD, tp); matHHC = myNmr.matConjugue(matHH); matC1 = myNmr.mult(matR2, matHH); matC2 = myNmr.mult(matHHC, matC1); matC3 = myNmr.mult(matC2, matRt); matC4 = myNmr.mult(matR, matC3); if (choix != 'A') amplitude[j] += (pop/coeff)*population[population.length - k - 1] *matC4.im[ligne][ligne - 1]; else switch (choiceSpin.getSelectedIndex()) { case 0: //spin 3/2 amplitude[j] += population[population.length - k - 1]*(2*matC4.im[2][1] + 2*Math.sqrt(3)*matC4.im[1][0])/coeff; break; case 1: //spin 5/2 amplitude[j] += population[population.length - k - 1]*(3*matC4.im[3][2] + 2*Math.sqrt(8)*matC4.im[2][1] + 2*Math.sqrt(5)*matC4.im[1][0])/coeff; break; case 2: //spin 7/2 amplitude[j] += population[population.length - k - 1]*(4*matC4.im[4][3] + 2*Math.sqrt(15)*matC4.im[3][2] + 2*Math.sqrt(12)*matC4.im[2][1] + 2*Math.sqrt(7)*matC4.im[1][0])/coeff; break; case 3: //spin 9/2 amplitude[j] += population[population.length - k - 1]*(5*matC4.im[5][4] + 2*Math.sqrt(24)*matC4.im[4][3] + 2*Math.sqrt(21)*matC4.im[3][2] + 2*Math.sqrt(16)*matC4.im[2][1] + 2*Math.sqrt(9)*matC4.im[1][0])/coeff; break; }//end of switch }//end of for j }//end of for k }//end of calculIntensite() private void afficheTransition() { choiceTransition.removeAll(); switch (choiceSpin.getSelectedIndex()) { case 0: //spin 3/2 choiceTransition.add("Central Transition"); choiceTransition.add("Satellite Transitions"); choiceTransition.add("All the Transitions"); break; case 1: //spin 5/2 choiceTransition.add("Central Transition"); choiceTransition.add("Inner Satellite Transitions"); choiceTransition.add("Outer Satellite Transitions"); choiceTransition.add("All the Transitions"); break; case 2: //spin 7/2 choiceTransition.add("Central Transition"); choiceTransition.add("Inner Satellite Transitions"); choiceTransition.add("Medium Satellite Transitions"); choiceTransition.add("Outer Satellite Transitions"); choiceTransition.add("All the Transitions"); break; case 3: //spin 9/2 choiceTransition.add("Central Transition"); choiceTransition.add("Inner Satellite Transitions"); choiceTransition.add("Medium Satellite Transitions"); choiceTransition.add("Next Satellite Transitions"); choiceTransition.add("Outer Satellite Transitions"); choiceTransition.add("All the Transitions"); break; }//end of switch choiceTransition.select(0); choix ='C'; }//end of afficheTransition() private void normalisationSignal() { switch (choiceSpin.getSelectedIndex()) { case 0: ComplexMat.setTailleMatrice (4); coeff = 5.0; switch (choix) { case 'C': ligne = 2; pop = Math.sqrt(4); break; case 'I': ligne = 1; pop = 2*Math.sqrt(3); break; } break; case 1: ComplexMat.setTailleMatrice (6); coeff = 35.0/2.0; switch (choix) { case 'C': ligne = 3; pop = Math.sqrt(9); break; case 'I': ligne = 2; pop = 2*Math.sqrt(8); break; case 'O': ligne = 1; pop = 2*Math.sqrt(5); break; } break; case 2: ComplexMat.setTailleMatrice (8); coeff = 42.0; switch (choix) { case 'C': ligne = 4; pop = Math.sqrt(16); break; case 'I': ligne = 3; pop = 2*Math.sqrt(15); break; case 'M': ligne = 2; pop = 2*Math.sqrt(12); break; case 'O': ligne = 1; pop = 2*Math.sqrt(7); break; } break; case 3: ComplexMat.setTailleMatrice (10); coeff = 165.0/2.0; switch (choix) { case 'C': ligne = 5; pop = Math.sqrt(25); break; case 'I': ligne = 4; pop = 2*Math.sqrt(24); break; case 'M': ligne = 3; pop = 2*Math.sqrt(21); break; case 'N': ligne = 2; pop = 2*Math.sqrt(16); break; case 'O': ligne = 1; pop = 2*Math.sqrt(9); break; } break; }//end of switch }//end of normalisationSignal() /**Start the applet*/ public void start() { //********* //buttonRun //********* buttonRun.addActionListener(new java.awt.event.ActionListener() { public void actionPerformed(ActionEvent e) { normalisationSignal(); /*transform pulse durations from microsecond unit to millisecond unit, because the interactions are in kHz unit*/ dureeFin = 0.001*Double.valueOf(textFieldDerniereDuree.getText()).doubleValue(); dureePremiere = 0.001*Double.valueOf(textFieldPremiereDuree.getText()).doubleValue(); dureePas = 0.001*Double.valueOf(textFieldPas.getText()).doubleValue(); nbou = (int)((dureeFin - dureePremiere)/dureePas); amplitude = new double[nbou + 1]; //FID intensity listExpIntensity = new double[nbou + 1]; //experimental intensity calculIntensite(Double.valueOf(textFieldOmegaQ.getText()).doubleValue(), Double.valueOf(textFieldOmegaRF.getText()).doubleValue()); scrollPaneSimulation.setScrollPosition(0, 0); //top part of scrollPane //*******show simulated line intensity versus the pulse duration******** double tp = 0; String s = "", s2 = ""; textAreaIntensite.setText(""); //clear intensity textAreaDuree.setText(""); //clear pulse duration for (int i = 0; i <= nbou; i++) { tp = dureePremiere + i*dureePas; s = textAreaIntensite.getText() + (int)(amplitude[i]) + "\n"; textAreaIntensite.setText(s); s2 = textAreaDuree.getText() + Math.round(100000*tp)/100.0 + "\n"; textAreaDuree.setText(s2); }//end of for i //********************************************************************** panelFit.setVisible(true); panelFitEstim.setVisible(true); } }); //********* //buttonFit //********* buttonFit.addActionListener(new java.awt.event.ActionListener() { public void actionPerformed(ActionEvent e) { /*thread for starting the simplex procedure with run() and for plotting the fitting curves*/ Thread th = null; if (buttonFit.getLabel() == "Fit") { //convert experimental intensities into numbers String s = textAreaExpIntensite.getText(), s1 = ""; int spaceAt = 0, startingForm = 0, ip = 0; while (true) { spaceAt = s.indexOf("\n", startingForm); if (spaceAt == -1) break; s1 = s.substring(startingForm, spaceAt); listExpIntensity[ip] = Double.valueOf(s1).doubleValue(); startingForm = spaceAt + 1; ip++; } if (ip != nbou + 1) { //check the number of data textAreaParameter.setText(""); textAreaParameter.append("The number of experimental" + "\n"); textAreaParameter.append("intensities is wrong!!!" + "\n"); textAreaParameter.append("May be type enter key" + "\n"); textAreaParameter.append("for the last intensity data!!!"); } else { //start fitting panelFit.setLocation(300, 0); setBackground(Color.pink); panelTransition.setVisible(false); panelSpin.setVisible(false); panelBlanc.setBounds(520 ,0, 240, hautPanel); add(panelBlanc); panelBlanc.setVisible(true); panelBlanc.setBackground(Color.white); buttonFit.setLabel("StopFit"); buttonPrevious.setEnabled(false); textAreaParameter.setText(""); double alpha = Double.valueOf(choiceSimplex.getSelectedItem()).doubleValue(); double val1 = Double.valueOf(textFieldOmegaQEst.getText()).doubleValue(); double val2 = Double.valueOf(textFieldOmegaRFEst.getText()).doubleValue(); //find the maximum values in listExpIntensity[] and in amplitude[] double maxExp = 0, maxCal = 0; for (int i = 0; i < nbou + 1; i++) if (Math.abs(listExpIntensity[i]) > maxExp) maxExp = Math.abs(listExpIntensity[i]); for (int i = 0; i < nbou + 1; i++) if (Math.abs(amplitude[i]) > maxCal) maxCal = Math.abs(amplitude[i]); double val3 = maxExp/maxCal; //scaling choiceSimplex.setEnabled(false); textFieldOmegaRFEst.setEnabled(false); textFieldOmegaQEst.setEnabled(false); textAreaExpIntensite.setEnabled(false); normalisationSignal(); //initial parameters for simplex procedure matEstimation[0][0] = val1; matEstimation[0][1] = val2; matEstimation[0][2] = val3; matEstimation[1][0] = alpha*val1; matEstimation[1][1] = val2; matEstimation[1][2] = val3; matEstimation[2][0] = val1; matEstimation[2][1] = alpha*val2; matEstimation[2][2] = val3; matEstimation[3][0] = val1; matEstimation[3][1] = val2; matEstimation[3][2] = alpha*val3; th = new Thread(OnePowderJDK118.this); th.start(); buttonGraph.setLabel("DeactivateGraph"); buttonGraph.setVisible(true); buttonGraph.setEnabled(true); labelAdvice.setVisible(true); showCurve = true; stopFit = false; } } else{ //stop fitting stopFit = true; th = null; buttonPrevious.setEnabled(true); buttonFit.setLabel("Fit"); choiceSimplex.setEnabled(true); textFieldOmegaRFEst.setEnabled(true); textFieldOmegaQEst.setEnabled(true); textAreaExpIntensite.setEnabled(true); buttonGraph.setEnabled(false); labelAdvice.setVisible(false); } } }); //************** //buttonPrevious //************** buttonPrevious.addActionListener(new java.awt.event.ActionListener() { public void actionPerformed(ActionEvent e) { panelSpin.setVisible(true); panelBlanc.setVisible(false); panelTransition.setLocation(220, 0); panelTransition.setVisible(true); panelFit.setLocation(440, 0); buttonPrevious.setEnabled(false); setBackground(Color.white); panelSpin.setBackground(Color.pink); panelFit.setBackground(Color.orange); panelTransition.setBackground(Color.pink); buttonGraph.setVisible(false); buttonFit.setLabel("Fit"); buttonRun.setEnabled(true); } }); //*********** //buttonGraph //*********** buttonGraph.addActionListener(new java.awt.event.ActionListener() { public void actionPerformed(ActionEvent e) { if (buttonGraph.getLabel() == "DeactivateGraph") { showCurve = false; buttonGraph.setLabel("ActivateGraph"); labelAdvice.setVisible(false); } else { showCurve = true; buttonGraph.setLabel("DeactivateGraph"); labelAdvice.setVisible(true); } } }); //********** //choiceEta //********** choiceEta.addItemListener(new java.awt.event.ItemListener() { public void itemStateChanged(ItemEvent e) { switch (choiceEta.getSelectedIndex()) { case 0: population = POPULATION_ETA00; break; case 1: population = POPULATION_ETA01; break; case 2: population = POPULATION_ETA02; break; case 3: population = POPULATION_ETA03; break; case 4: population = POPULATION_ETA04; break; case 5: population = POPULATION_ETA05; break; case 6: population = POPULATION_ETA06; break; case 7: population = POPULATION_ETA07; break; case 8: population = POPULATION_ETA08; break; case 9: population = POPULATION_ETA09; break; case 10: population = POPULATION_ETA10; break; }//end of switch } }); //********** //choiceSpin //********** choiceSpin.addItemListener(new java.awt.event.ItemListener() { public void itemStateChanged(ItemEvent e) { textAreaIntensite.setText(""); //clear intensity textAreaDuree.setText(""); //clear pulse duration afficheTransition(); } }); //**************** //choiceTransition //**************** choiceTransition.addItemListener(new java.awt.event.ItemListener() { public void itemStateChanged(ItemEvent e) { String rangTransition = choiceTransition.getSelectedItem(); if (rangTransition == "Central Transition") choix = 'C'; else if (rangTransition == "Satellite Transitions") choix = 'I'; else if (rangTransition == "Inner Satellite Transitions") choix = 'I'; else if (rangTransition == "Medium Satellite Transitions") choix = 'M'; else if (rangTransition == "Next Satellite Transitions") choix = 'N'; else if (rangTransition == "Outer Satellite Transitions") choix = 'O'; else if (rangTransition == "All the Transitions") choix = 'A'; textAreaIntensite.setText(""); //clear intensity textAreaDuree.setText(""); //clear pulse duration } }); }//end of start() /**Stop the applet*/ public void stop() { } /**Destroy the applet*/ public void destroy() { } /**Get Applet information*/ public String getAppletInfo() { return "Applet Information"; } /**Get parameter info*/ public String[][] getParameterInfo() { return null; } private double mySimplexFonction(double omegaQFit, double omegaRFFit, double norme) { double sqsum = 0; calculIntensite(omegaQFit, omegaRFFit); for (int i = 0; i <= nbou; i++) { double tmp = norme*amplitude[i] - listExpIntensity[i]; sqsum += tmp*tmp; } return sqsum; }//end of mySimplexFonction() private void plotCourbe(double b0, double b1, double b2) { //plot experimental intensity data and the fitting curve double inc; Graph("One pulse to powder", "tp(µs)", "S(tp)"); setSymbol(true); setLine(false); setColor(blue); setTitle("experimental"); for (int i = 0; i <= nbou; i++) { inc = Math.round(100000*(dureePremiere + i*dureePas))/100.0; add(inc, listExpIntensity[i]); } nextGraph(); setColor(black); setTitle("fit"); calculIntensite(b0, b1); for (int i = 0; i <= nbou; i++) { inc = Math.round(100000*(dureePremiere + i*dureePas))/100.0; add(inc, b2*amplitude[i]); } showGraph(); }//end of plotCourbe() public void run() { //run the simplex procedure /* The flexible simplex method allows you to find the minimum for a function * with multiple variables. The method obtains the minimum for a function * with N variables by examining the function values at N + 1 points. * This method locates the points with the best and worst values * then attempts to replace the worst-value point with a better point. * This replacement process involves expanding and contracting the simplex * near the worst-value point to determine a better replacement point. * The iterations of the method shrinks the multidimensional simplex around * the minimum point. * Given * # N variables that evaluate function f(X) * # The matrix X, which contains N + 1 rows and N columns * # The array of function values, Y(containing N + 1 values) * # The maximum number of iterations, MAXITER * # The reflection factor, which is a positive value that may be unity * # The expansion factor, which is greater than one (usually two) * # The contraction factor, a value between zero and one (ususally 0.5) */ //(1) set Iter = 0 //(2) set P = N + 1 //the 4 constants used in simplex procedure final int REFLECTFACT = 1; final int EXPANDFACT = 2; final double CONTRACTFACT = 0.5; final int MAXITER = 1000; int numIter = 0, numPoints = numVars + 1, worstI = 0, bestI = 0; double y1 = 0, y2 = 0, x0 = 0; //(5) set the convergence flag to false boolean flag = false; double y[] = new double[numVars + 1]; //(3) allocate dynamic arrays x1, x2 and centroid to contain N elements double x1[] = new double[numVars]; double x2[] = new double[numVars]; double centroid[] = new double[numVars]; //(4) calculate the values for y() using x1() for (int i = 0; i < numPoints; i++) { for (int j = 0; j < numVars; j++) x1[j] = matEstimation[i][j]; y[i] = mySimplexFonction(x1[0], x1[1], x1[2]); } //************************* boucle while ********************************* //(6) repeat the next step while Iter < M and the convergence flag is false while ((numIter < MAXITER) && (!stopFit)) { String sx = ""; textAreaParameter.setText(""); sx = "Iteration: " + numIter + " \n"; textAreaParameter.append(sx); numIter += 1; //(6.1) increment Iter //(6.2) find the indices of the worst A and best points worstI = 0; bestI = 0; for (int i = 1; i < numPoints; i++) { if (y[i] < y[bestI]) bestI = i; else if (y[i] > y[worstI]) worstI = i; } //provide the fitting parameters sx = Math.round(1000*matEstimation[bestI][1])/1000.0 + " ( " + Math.round(1000*matEstimation[worstI][1])/1000.0 + " ) kHz\n"; textAreaParameter.append("RF field: " + sx); sx = Math.round(1000*matEstimation[bestI][0])/1000.0 + " ( " + Math.round(1000*matEstimation[worstI][0])/1000.0 + " ) kHz\n"; textAreaParameter.append("QCC: " + sx); sx = Math.round(1000*matEstimation[bestI][2])/1000.0 + " ( " + Math.round(1000*matEstimation[worstI][2])/1000.0 + " )\n"; textAreaParameter.append("Scale: " + sx); sx = y[bestI] + "\n"; textAreaParameter.append(sx); sx = y[worstI] + "\n"; textAreaParameter.append(sx); if (showCurve) { try{Thread.sleep(10);} catch (InterruptedException e) {} //plotCourbe should be located after try otherwise problem with Netscape plotCourbe(matEstimation[worstI][0], matEstimation[worstI][1], matEstimation[worstI][2]); } //(6.3) calculate centroid (exclude the worst point A) //= center of gravity G for (int i = 0; i < numVars; i++) { centroid[i] = 0; for (int j = 0; j < numPoints; j++) { if (j != worstI) centroid[i] += matEstimation[j][i]; } centroid[i] /= numVars; } //(6.4) calculate reflected point A' //(symmetrical point of A relative to G): A' = 2*G - A for (int i = 0; i < numVars; i++) x1[i] = (1 + REFLECTFACT) * centroid[i] - REFLECTFACT * matEstimation[worstI][i]; //(6.5) set y1 = function value at point x1 or f(A') y1 = mySimplexFonction(x1[0], x1[1], x1[2]); //************************** beginning of if **************************** //(6.6) if y1 or f(A') < best y, that is, //f(A') is the best value found so far if (y1 < y[bestI]) { //perform the following tasks: //(6.6.1) calculate the expandeded point x2 or A" = 2A' - G: //a point twice as far from G as A' for (int k = 0; k < numVars; k++) x2[k] = (1 + EXPANDFACT) * x1[k] - EXPANDFACT * centroid[k]; //(6.6.2) set y2 = function value at point x2, or f(A") y2 = mySimplexFonction(x2[0], x2[1], x2[2]); //(6.6.3) if y2 or f(A") < best y if (y2 < y[bestI]) { //then replace worst point matEstimation //with point x2 or A" and resume at step 7 for (int m = 0; m < numVars; m++) matEstimation[worstI][m] = x2[m]; } else { //else replace worst point matEstimation with point x1 or A' //and resume at step 7 for (int n = 0; n < numVars; n++) matEstimation[worstI][n] = x1[n]; } } //******************************* else *********************************** else { //(6.7) y1 or f(A') >= best y flag = true; for (int q = 0; q < numPoints; q++) { if ((q != worstI) && (y1 <= y[q])) { flag = false; //f(A') is between best y and the others y, except A break; //exit for } } if (flag) {//y1 or f(A') is not between best y and the others y, except A if (y1 < y[worstI]) { //f(A') is between best y and f(A), //(6.7.1) replace worst point by x1 or A' for (int p = 0; p < numVars; p++) matEstimation[worstI][p] = x1[p]; y[worstI] = y1; } //f(A') >= f(A) //(6.7.2) calculate contracted point x2 //or B = (A + G)/2 or B = (A' + G)/2 for (int s = 0; s < numVars; s ++) x2[s] = CONTRACTFACT * matEstimation[worstI][s] + (1 - CONTRACTFACT) * centroid[s]; //(6.7.3) set y2 = function value at point x2 or f(B) y2 = mySimplexFonction(x2[0], x2[1], x2[2]); if (y2 > y[worstI]) { //f(B)>f(A) or f(B)>f(A'), //ALL WE DID ARE NOT USEFUL //store best matEstimation in x2 for (int u = 0; u < numVars; u++) x2[u] = matEstimation[bestI][u]; //and build a new simplex obtained by dividing //all edges leading to the point yielding the best value by 2 for (int j = 0; j < numPoints; j++) for (int q = 0; q < numVars; q++) matEstimation[j][q] = 0.5 * (x2[q] + matEstimation[j][q]); } //resume 7 else { //(else of y2), f(B)