在 由离散点绘制的曲线上 找出零点对应的横坐标

chenjc18

针对ChaChing 兄 对上述程序提出的见解
我将问题暂时解决了一下：
请看程序：
此程序只是使用于含有两个零(及以上)的离散序列，对于含有1个或者连续出现两个零点的情形 可以按照我提出的 分段求解 思路 对下面的程序稍加修改 即可使用。

请大家，试试看，还有什么问题。
鉴于我刚刚学习MATLAB，其中出现幼稚的语句也在所难免，还请大家多提出建议。谢谢了

1. clear;close;
   
2. 
   
3. y=[1 2 -3 -1 0 2 8 -1 -1 4 0 -1 1 -4 7 -2];
   
4. x=linspace(0,20,length(y));
   
5. zero=zeros(length(x),1);
   
6. x=x';
   
7. y=y';
   
8. plot(x,y,'b',x,zero,'r');
   
9. ind=find(y==0);
   
10. X1=x(ind); %找到自身就有的零点对应的x
    
11. %段一 从段首到第一个自身的零点
    
12. %以下是对段一求解零点
    
13. a=y(1:(ind(1)-1));
    
14. b=find(a(1:end-1).*a(2:end)<0);
    
15. c=[a(b) a(b+1)];
    
16. d=[x(b) x(b+1)];
    
17. for k=1:length(b)
    
18. X2(k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
19. end
    
20. %中间的各段
    
21. clear a b c d ;
    
22. for i=1:(length(ind)-1)
    
23. if((ind(i+1)-ind(i)~=1)|(ind(i+1)-ind(i)~=2))
    
24. a=y((ind(i)+1):(ind(i+1)-1));
    
25. a1=x((ind(i)+1):(ind(i+1)-1));
    
26. b=find(a(1:end-1).*a(2:end)<0);
    
27. c=[a(b) a(b+1)];
    
28. d=[a1(b) a1(b+1)];
    
29. for k=1:length(b)
    
30. X_3(i,k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
31. end
    
32. else continue;
    
33. end
    
34. end
    
35. X_30=reshape(X_3,1,[]); %矩阵转向量；
    
36. X3=X_30(find(X_30)); %剔除系统自动补的零
    
37. %最后一段，最后一个自身的零点 到 y序列末
    
38. clear a b c d a1;
    
39. a=y((ind(length(ind))+1):length(y));
    
40. a1=x((ind(length(ind))+1):length(y));
    
41. b=find(a(1:end-1).*a(2:end)<0);
    
42. c=[a(b) a(b+1)];
    
43. d=[a1(b) a1(b+1)];
    
44. X=zeros(size(b));
    
45. for k=1:length(b)
    
46. X4(k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
47. end
    
48. X_1=cat(2,X1',X2);
    
49. X_2=cat(2,X3,X4);
    
50. X=cat(2,X_1,X_2);
    
51. X=sort(X) %按升序排列 得到零点序列
    

复制代码

[ 本帖最后由 chenjc18 于 2009-1-5 21:36 编辑 ]

ChaChing

楼主了解涵义了! 人不适, 明天再看看!

chenjc18

我将16楼的做了改到 可以处理任何形式的离散序列：适用于：不含零点 含有一个零点 含有两个（不）连续零点 含有两个以上零点
请大家，看一下程序，有什么问题
谢谢了
我写的这段程序可能有点太乱了，嘿嘿，我已经尽力 欢迎大家批评指正，大家共同进步！

1. 
   
2. clear;close;
   
3. y=[1 2 -3  -1  0 1 -1  0 2 0 8 0 -1 -1  4  1 -4  7 -2];
   
4. x=linspace(0,20,length(y));
   
5. zero=zeros(length(x),1);
   
6. x=x';
   
7. y=y';
   
8. plot(x,y,'b',x,zero,'r');
   
9. ind=find(y==0); %找出自身零点所在下标
   
10. if(isempty(ind)) %如果原离散序列不含有零，则进行程序一
    
11. %%%%%%%%程序一%%%%%%%%%
    
12.     b=find(y(1:end-1).*y(2:end)<0);
    
13.     c=[y(b) y(b+1)];
    
14.     d=[x(b) x(b+1)];
    
15.     for k=1:length(b)
    
16.         X(k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
17.     end
    
18.     X %输出x
    
19. else  %否则执行程序二
    
20. %%%%%%%%%程序二%%%%%%%%%%%  
    
21. X1=x(ind);   %找到自身就有的零点对应的x
    
22. %%段一 从段首到第一个自身的零点
    
23. a=y(1:(ind(1)-1));
    
24. b=find(a(1:end-1).*a(2:end)<0);
    
25. c=[a(b) a(b+1)];
    
26. d=[x(b) x(b+1)];
    
27. X=zeros(size(b));
    
28. for k=1:length(b)
    
29. X2(k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
30. end
    
31. %中间的各段
    
32. clear a b c d ;
    
33. if((length(ind)==2)&(ind(2)-ind(1)==1)|(length(ind)==1))
    
34.     X_3=[0];%如果原离散序列含有一个，或者连续两个零点，则执行此句
    
35. else        %如果含有两个（不连续）及以上个零点 执行下面的语句
    
36. for i=1:(length(ind)-1)
    
37.     if((ind(i+1)-ind(i)~=1)|(ind(i+1)-ind(i)~=2))%如果出现相邻的或间隔为1的零点，则跳过
    
38.         a=y((ind(i)+1):(ind(i+1)-1));            
    
39.         a1=x((ind(i)+1):(ind(i+1)-1));
    
40.         b=find(a(1:end-1).*a(2:end)<0);
    
41.         c=[a(b) a(b+1)];
    
42.         d=[a1(b) a1(b+1)];
    
43.         for k=1:length(b)
    
44.         X_3(i,k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
45.         end
    
46.     else continue;
    
47.     end
    
48. end
    
49. end
    
50. X_30=reshape(X_3,1,[]); %矩阵转向量；
    
51. X3=X_30(find(X_30));    %剔除系统自行添加的零
    
52. %%最后一段，最后一个自身的零点 到 y序列末
    
53. clear a b c d a1;
    
54. a=y((ind(length(ind))+1):length(y));
    
55. a1=x((ind(length(ind))+1):length(y));
    
56. b=find(a(1:end-1).*a(2:end)<0);
    
57. c=[a(b) a(b+1)];
    
58. d=[a1(b) a1(b+1)];
    
59. for k=1:length(b)
    
60. X4(k)=interp1([c(k,1), c(k,2)],[d(k,1),d(k,2)],0);
    
61. end
    
62. %%将各段计算得出零点放在一起
    
63. X_1=cat(2,X1',X2);      
    
64. X_2=cat(2,X3,X4);
    
65. X=cat(2,X_1,X_2);
    
66. X=sort(X)                %按升序排列
    
67. end
    

复制代码

ch_j1985

LZ还真有研究的劲头，佩服！

ChaChing

真服了楼主的态度! 欣赏!
还记得7F的回覆吗!? 其中说过interp1好像根本不需要!
18F的程序太长了, 逻辑思绪个人感觉太复杂了, 有点懒得看, 不好意思!
欣赏楼主, 列上我的, 不是最佳, 参考参考!
clc; clear;
y=[1 2 -3  -1  0 1 -1  0 2 0 8 0 -1 -1  4  1 -4  7 -2];
x=linspace(0,20,length(y)); x=x'; y=y'; xx=[]; plot(x,y); grid on;
for k=1:length(y)-1
   if y(k)==0, xx=[xx; x(k)];
   elseif y(k)*y(k+1)<0, xx=[xx; x(k)-y(k)*(x(k+1)-x(k))/(y(k+1)-y(k))]; end
end
if y(end)==0, xx=[xx; x(end)]; end
xx

其中x(k)-y(k)*(x(k+1)-x(k))/(y(k+1)-y(k))亦可用interp1([y(k),y(k+1)],[x(k),x(k+1)],0)

[ 本帖最后由 ChaChing 于 2009-1-6 16:09 编辑 ]

chenjc18

谢谢 ChaChing
您编写的程序真棒！
又长了很多知识，谢谢了